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SECTION  I 
INTOODUCnCN 

Accurate  prediction  of  free-f light  aerodynandc 
characteristics  is  of  the  utmost  inportance  to  the  weapon 
airframe  designer.  Mission  effectiveness  is  highly  dependent  on 
his  ability  to  design  an  airframe  with  desirable  aerodynamic 
qualities  to  deliver  the  payload  on  target  with  optimum 
maneuvering  performance.  Often  this  is  a  time-consuming 
iterative  effort  coupling  various  aerodynamic  predictive  methods 
and  ejgjensive  wind  tunnel  testing  for  airframe  validation. 

As  a  result,  the  Air  Force  Armament  Laboratory  has  embarked 
on  a  multi-year  effort  to  derive  accurate  computational  fluid 
dynamic  techniques  to  enhance  weapon  airframe  design  procedures 
and  reduce  wind  tunnel  validation  testing.  The  first  product 
from  this  effort  is  Program  EACILE  -  Flow  Solver  (MISSE);  a  multi¬ 
block,  implicit,  steady-state  Euler  algorithm.  Program  EAGLE  has 
been  designed  to  aid  in  the  analysis  of  both  freestream  and 
interference  aerodynamic  characteristics  of  advanced  arbitrarily 
shaped,  finned  or  unfinned  configurations  in  the  catpressible 
subsonic,  transonic  and  supersonic  Mach  ranges  at  relatively  low 
incidence  angles. 

This  technical  report  is  designed  to  be  a  user's  manual  for 
the  Program  EAGLE  -  Flow  Solver.  In  this  report  an  indepth 
discussion  is  given  on  the  numerical  solution  of  the  three- 
dimensional  Euler  equations  (Section  II  )  .  In  addition,  a  section 
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r 


is  provided  that  details  the  entire  set  of  input  required  by  the 
Flow  Solver  including  all  default  values  (Section  III).  A 
detailed  discussion  of  the  main  routine  (MISSE),  each  of  the  39 
subroutines  and  the  interaction  between  all  portions  of  the  code 
(Section  IV),  is  also  included  for  the  user.  And  finally,  in 
Section  V  (Applications)  a  ccnplete  discussion  of  several 
aerodynamic  applications  is  given  to  familiarize  the  user  with 
the  detailed  workings  of  the  Program  EAGLE  -  Flow  Solver. 
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SECTION  II 
IHEDREnCAL  ASPECTS 

1 .  GOVn^ING  EQUATIONS 


The  three-dimensional,  time  dependent  Euler  equations 
written  in  strong  conservation  law  form  are 


(a^/a't)  + 

(aVay)  +  (ah/a'^)  = 

0 

where 

<■'  1 

p 

P'  ~ 

V 

U'  ! 

i  ) 

ipc> 

i 

p^i^ 

■  pCi^r 

q  =;  |iv  !  ,  1  = 

A  A 

puv 

,  g  =  pv^+P 

A 

,  h 

P s  /\  A 

=  .  pvw 

puw 

A  A  /•, 

pvw 

.  A  2  A  , 

pw  +P  . 

:  A 

e 

; _ 1 

u(e+P) 

and 

^  =  (P/y-1)  +  (p/2)[^iW+^^] 

Here  p  is  the  mass  density  of  the  ocitpressible  fluid,  li,  v, 
cind  w  are  the  velocity  corrponents  in  the  x,  y,  and  z  Cartesian 
coordinate  directions,  respectively,  P  is  the  pressure,  and  e  the 
total  specific  energy  of  the  fluid.  The  relationship  between  e 

and  P  is  the  result  of  making  the  perfect  gas  assunption,  and  y 
is  the  ratio  of  specific  heats.  Carets  on  the  quantities  in  Equation 
( 1 )  are  used  to  indicate  that  dimensional  quantities  are  being 
used.  One  way  to  non-dimensionalize  Equation  (1 )  is  to  define  the 
following  non-dimensional  quantities: 

X  =  x/lr  '  y  =  Wir  '  z  =  z/lj.  ,  t  =(a,,-,t  )/l  j. 

u  =  u/a^^,,  V  =  v/a^,  w  = 

p  =  |i/p^^^,  e  =  e/(,,^^  a2,)  ,  P  =  ) 

1  /2 

where  I  is  any  convenient  reference  length,  a=  (yP,.,/|i  ^  '  the 

'  (X> 
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speed  of  sound  in  undisturbed  gas,  and  n  is  the  density  of 

^00 

undisturbed  gas.  After  some  manipulation,  the  non-dimensional 
Euler  equations  may  be  written  in  exac±ly  the  same  form  as  Equation 
( 1 )  with  all  carets  dropped. 

lb  simplify  numerical  treatment  of  boundary  conditions  around 
general  geometries  the  Euler  equations  are  recast  in  terms  of 
general  boundary  conforming  curvilineaur  coordinates.  The 
curvilinear  coordinates  are  defined  as 


^  =  t-(x,y,2) 

7^  =  |](x,y,z) 

T  =  t 

Applying  the  transformation  Equation  (3)  to  the  non-dimensional 
Euler  equations  (Equation  (1 )  with  all  carets  dropped)  gives 

(.30/31)  +  (3F/3^')  +  (3G/3ii)  +  (.3H/3C)  =  o 

where 


Q  =  J 


P 

pu 
pv 
pw 
e  , 


pu 

puU  + 


F  =  J  pvU  + 

pwU  +  ^^P 
|j[J(e+P’) 

1  pv 

puV  +  ■p^P 

G  =  J  ;  pvv  +  p^p 


and  U  =  tu  +  tv  +  )-’w 
Sx  '»y  '»z 


pwV  +  p  P 
V(efP) 


and  V 


Hv''  V  +  p  w 


H  =  J 


pW 

puW  +  ^^P 
pww  +  c  P 


and  W  -  u  +  r  V  +  r  w 
•  X  •  y  •  z 


_W(efP) 

in  the  above  equation,  J  is  the  Jacobian  of  the  inverse 


(3) 


(4) 
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transforJtaticyi.  i-e.  3(x,y,z)/9(t-,^|,^ )  and  is  given  by 
The  metric  quantities  are  given  by 


Sy  ~  ^  11y  ~  (x^-Z^-2^X^  ),  =  J  ^  (z^Xj^-X|;Z^^  ) 


The  details  of  this  transformation  are  given  in  Reference  1 • 


2.  FLUX-VBCIDR  SPLITTING  ( Steger-Warming ) 

Hyperbolic  partial  differential  equations,  such  cis  the  Euler 
equations,  are  characterized  by  the  existence  of  a  limited  domain 
of  dependence.  Tlie  solution  at  a  point  does  not  depend  on  every 
other  point  in  the  field;  this  means  that  information  travels 
only  in  certain  characteristic  directions.  Numerical  schemes 
intended  to  solve  hyperbolic  equations  are  usually  enhaix:»i  by 
insuring  that  the  numerical  method  prqpagates  information  in  the 
direction  specified  by  the  partial  differential  equation.  This 
can  be  done  by  using  an  upwind  method,  or  one  in  which  the 
difference  operator  is  taken  in  the  direction  from  vrtiich  the 
information  should  ccme.  Stability  properties  are  often  improved 
by  upwinding,  and  it  is  usually  unnecessary  to  add  smoothing 
terms  or  artificial  viscosity  to  an  upwind  method. 

The  three-dimensional  Euler  equations,  liquation  (4),  are  a 
hyperbolic  system  of  five  equations  and  hence  have  five 
characteristic  velocities  in  each  of  the  three  spatial 
directions.  These  characteristic  velocities  are  determined  from 
the  quasilinear  form  of  Equation  ( 4 ) , 
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(QQ/ai)  +  A(3Q/a^)  +  B(,3Q/3li)  +  c(3Q/ap  =  0 


(5) 


where  the  matrices  A,  B,  and  C  are  given  by 

A=  (3F/3Q)  ,  B=  (3G/3Q)  ,  C=  (3H/3Q)  .  (6) 

The  eigenvalues  of  A  are  the  characteristic  velocities  in  the  J;' 
direction,  the  eigenvalues  of  B  are  the  characteristic  velocities 
in  the  direction,  and  similarly,  the  eigenvalues  of  C  are  the 
characteristic  velocities  in  the  ^  direction. 

Since  F,  G,  and  H  are  identical  except  that  where  ^  appears 
in  F,  will  afpear  in  G  and  ^  will  appeau:  in  H,  extra  writing 
can  be  avoided  by  letting  K  represent  either  F,  G,  or  H,  when  k 
represents  either  ^  or  ^  respectively.  Then  define 

K  =  {3K/3Q)  (7) 

which  corresponds  to  A,  B,  or  C  depending  on  the  meaning  of  K. 

The  general  flux  Jacobian  matrix  K  and  its  eigenvalues  are 
derived  in  Reference  2.  The  eigenvalues  of  the  matrix  K  are 

~  ^k  ~  ^k  ~  ^y'^  ^  ~  0Jq 

~  0k  ^  ”  0]^  "-^^x  ^y  ^  ^ 

X^  =  0k  ~ 

where  c  is  the  speed  of  sound. 

It  is  possible  to  split  the  flux  vector  K  into  three  parts, 
one  corresponding  to  each  of  the  distinct  eigenvalues  of  K  given 
above  (for  the  details  of  this  splitting,  see  Reference  2). 

The  flux  vector  K  is  then  written  as 

K  =  Xk^i  + 


where 
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K 


1 


1 

J(y-1)/y'  pv 
,  pw 


(p/2)(u2 


+ 


=  J/(2y  ) 


Ipu  +  pd^  ; 
;  pv  +  pcky 

pw  +  pck2  ^  ' 
|_e  +  P  + 


(10a) 


(10b) 


I 

I  pu  -  pckx  I 

K5  =  J/(2y  )  !  pv  -  pcky  ; 

pw  -  pdc2 
ie  +  P  -  pc0,^; 

and 

kjf  =  (kx/|yk|  )  =  kx/(kx  +  ky  +  k2)’/^ 
ky  =  (ky/|yk|  ) 
k2  =  (k2/|7ki  ) 

8j^  =  kjjU  +  iCyV  +  jc2W. 

Ibe  sign  of  in  Equation  (8)  determines  from  which  direction 
information  should  be  used  to  determine  the  corresponding  portion 
of  flux,  K|  for  1  =1,4,  and  5. 

The  discussion  above  provides  a  rationale  for  writing  the 
flux  vector  K  eis  the  sum  of  two  vectors,  and  K~; 

K  =  +  K~ 

K"*"  is  associated  with  the  eigenvalues  that  have  positive  signs 
and  K~  is  associated  with  the  eigenvalues  that  have  negative 

signs.  The  specific  method  used  to  calculate  K'*'  and  K“  will  be 
discussed  in  the  next  section. 


(10c) 


(lOd) 
(lOe) 
(lOf ) 
(lOg) 


(11) 
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a.  Finite  Volume  Discaretization 
A  finite  volume  discretization  of  Equation  (4)  balances  the 
increase  of  the  conserved  quantity  in  a  ccnputational  cell,  or 
volume,  with  the  flux  of  the  quantity  through  the  surface  of  the 
cell .  Figure  ( 1 )  d^icts  a  portion  of  the  ccnputational  dcinain 
with  a  typical  cell  labeled.  Assume  the  dependent  variables  are 
constant  in  the  interior  of  cell  (i,  j,  k),  and  that  the  f.  « 
vectors  F,  G,  and  H  are  constant  over  the  constant  ^  / 

and  ^  surfaces  of  the  cell,  respectively.  Then,  an 
iitplicit  discretization  of  Bguation  (4)  is 

-  Q"i,j,ic'A5At)ar 

+  'G^'i,j.l/2,k  -°"^’i,3-1/2,I:>^5ACAT 

<>i"*'l,j,k.l/2  -  “"*'i,3,k-l/2’'^5Al|AT.  0 
This  can  be  written  as 


S  M  V 


(13) 


where 


for  exanple,  is  defined  by 


'^^^i,j,k  ~  ^  ^f-^i+1/2,  j,k  ^i-1/2,j,k^' 


cind  6^^  and  are  defined  analogDusly.  Letting  F"*",  F  ,  G'*’,  G  , 

and  H~  be  the  split  flux  vectors  for  F,  G,  and  H  as  given  by 
Equation  (11 ),  and  the  iirplicit  split  flux  discretization  is  given  by 


A(f  +  AtCS^CF"^  +  F")"+^  +  6^^(G^  +  =  0 


{ 14) 


The  fluxes  cure  nonlinear  functions  of  the  dependent 
variables,  and  must  be  linearized  to  obtain  a  system  of  equations 

that  can  easily  t5e  solved.  Ream  ani  Warmiivi'  .ml  Rriley  and 

McDonald^  have  used  the  linearization 
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kH+I  ^  j^n  ^  k"aQ"  +  O(At^) 


(15) 


wheare  =  (HK/HQ)*^  for  a  typical  flux  term  K.  The  flux  Jacobian 
matrices  required  to  linearize  Equation  (14)  are 

=  (aFVaQ)"^ ,  A“  =  (aF'/aQ)" , 

=  (acVao)"  /  B'  =  (aG“/aQ)” »  <’6) 

c+  =  (aH-^/aQ)" »  c"  =  (an'/aQ)" 

and  are  derived  in  Reference  1 •  Using  this  linearization,  Equation 
(14)  becomes 


[I  +  At(5\a^ 
s 


^  ^  V"  ^  'V'  ^ '  V"  ^ '  V' 


+  )]aq" 


(17) 


where  a  distinction  has  been  made  between  the  ijiplicit  spatial 
difference  operators  and  the  explicit  spatial  difference 
operators  by  using  superscripts  i  and  e  respectively.  The  dots 
are  used  to  indicate  that  the  differaice  operators  apply  to  the 

product  of  the  Jacobian  matrices  with  AQ^*  If  6  c' 


r  ^  11 


r 


are  chosen  such  that 


ft 


e 

11 


and  5 


i 

C' 

s 


b 


b 


i 
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=  (3/a^)  +  0(A^^), 

(18a) 

=  1<3/0ti  ^ 

(18b) 

=  (3/aC)  +  O(AC^), 

(18c) 

and  5^^  211^6  chosen  such  that 

=  id/'dV  ^ 

( 18d) 

=  (a/aii  ^  '^^A'd  1' 

( 18e) 

=  (6/aci  +  o(Aci. 

( 18f ) 

then  substitution  into  Equation  (17)  yields  a  method  that  is  second 
order  accurate  in  space  if  Aj;»  Al)'  ^nd  AC  ^re  0(j;i7  )  and 
first-order  accurate  in  time  when  used  with  Equation  (17).  The 
method  has  second  order  spatial  accuracy  when  converged  to 
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steady-state  since  approaches  0- 

For  this  algorithm,  the  inplicit  spatial  differences  need 
only  be  first  order  accurate  and  are  evaluated  using  a  one-point 
dependent  variable  extrapolation;  however,  the  esqjlicit 
differences  in  Egs.  (18)  <±ttain  second  order  accuracy  by 
evaluating  the  flux  at  cell  faces  from  two-point  d^3endent 
vcuriable  extrapolations. 

F'*’  is  evaluated  using  a  two-point  d^aendent  variable 
extrapolation  from  the  "positive”  side  and  a  two-point 
d^5endent  variable  extrapolation  from  the  "negative"  side  for  F”. 

For  a  particulcur  cell  face,  the  scheme  uses 

Qj+1/2  =  (3/2)Qj  -  (1/2)Qj_,  (19a) 

to  calculate  a  set  of  "positive"  eigenvalues,  (Q^),  and  left 
split  fluxes,  F|  (Q"*"),  and  uses 

Qj+1/2  =  <3/2)Qj^^  -  (l/2)Qj^2 

bo  calculate  a  set  of  "negative"  eigenvalues,  \*  (Q"),  and  right 
split  fltutes,  F”(Q~).  The  flux  at  the  cell  face  is  then  set  to 

=  Z  (1/2)[X'(Q‘^)  +  F,  (Q*^)  (20a) 

T 

and 

F‘  =  2  (1/2)[X'(Q~)  -  IX‘(Q")|]  F  (Q-).  (20b) 

I 

Note  that  if  a  "positive"  and  "negative"  eigenvalue  have 
different  signs,  then  the  corresponding  split  flux  may  either  be 
sunmed  from  both  sides,  if  the  sign  changes  from  positive  to 
negative,  or  from  neither  side,  if  the  sign  changes  from  negative 
to  positive. 

This  technique  seems  to  require  less  storage  and  executes 
more  quic)cly  than  previous  methods.  It  also  totally  eliminates 

1 1 


upstream  prcpagation  of  information  in  supersonic  regions. 
E>q)erience  with  this  method  has  shown  that  the  dispersive  effects 
cfcwnstream  of  shocks,  characteristic  of  previous  methods,  is  also 
reduced. 

b.  Aiproximate  Factorization 
Hie  finite  volume  schane  described  ty  Equation  (17)  is  not 
practical  to  use  because  the  system  of  algebraic  equations  that 
are  generated  have  a  very  large  bandwidth.  There  are  many 
possible  approximate  factorizations  of  this  method  that  can  be 
used  to  make  the  solution  process  easier,  as  noted  by  Steger  and 
Warming^,  for  exairple.  Several  factorizations  have  been 
recently  applied  by  Whitfield^  and  by  Anderson,  Thcmas,  and 

Whitfield  .  The  approximate  factorization  used  here  is 
given  by 

[I  +  +  aV.  +  6^c+0] 

[I  +  +  6V.  +  )]aq"  =  -AtR" 

where  r"  =  5®  e"  +  5®  g"  +  6®  h". 

^  11  C 

This  two-factor  scheme  is  very  attractive  in  that  it  consists  of 
the  solution  of  a  sparse  block  lower  triangular  system  followed 
by  the  solution  of  a  sparse  block  upper  triangular  system  of 
eq’jations  given  by 

[I  +  At<6Ja+.  +  x’  =  -AtR" 

5  11  C 

[I  +  At^S^a’  +  6^5“  +  hjc- ]  aq"  = 

?  11  T 

Solution  of  (21a)  is  done  by  a  simple  forward  substitution 
and  solution  of  (21b)  by  a  simple  back  SLibstitution.  A 
stability  analysis  for  several  factorization  schemes  was 
presented  by  Anderson,  et  al.,  in  Reference  8.  This  two-factor 


(21a) 

(21b) 


12 


scheme  was  found  to  be  least  sensitive  to  Courant  number  and  was 
stable  for  all  Courant  numbers  up  to  35,  which  was  the  maximum 
investigated. 

3.  FLUX-DIFFEKENCE  SPLITTING  (Roe)^ 

The  flux- vector  split  scheme  due  to  Steger  and  Warming^ 
described  in  Section  2.2  has  been  shown  to  capture  shocks  in  as 
few  as  two  mesh  points.  However,  a  flux-difference  split  scheme 

due  to  Roe^®  hcis  been  shown  to  capture  shocks  in  as 

few  as  zero  mesh  points.  Bram  van  Leer,  et.al.^^ 
recently  demonstrated  that  only  one-fourth  to  one-half  as  many 
points  were  required  to  resolve  a  shear  layer  using  a  first  order 
Roe  scheme  as  compared  to  the  number  of  points  required  by  a 
second  order  flux-vector  split  scheme.  To  allow  extension  of  the 
current  capabilities  of  the  Flow  Solver  to  a  partial  or  full 
Navier-Stokes  code,  a  flux-difference  split  scheme  has  been 
inplemented. 

a.  Approximate  Riemann  Solver  and  Roe  Averaging 
The  one-dimensional  Ebler  equations  are  considered  for 
sinplicity  since  in  the  one-dimensional  analysis  the  waves  move 
normal  to  the  cell  interface.  In  a  finite-volume  formulation  for 
multi-dimensional  flow  the  assumption  is  made  that  all  waves  move 
normal  to  the  interfaces  and  the  approximate  Riemann  solver  is 
applied  to  each  interface  based  on  this  assumption. 

The  one-dimensional  EXiler  equations  are: 
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(ao/at)  +  (aF/ax)  =  o 


(22) 


Roe's  idea  was  to  approximate  the  Riamann  problem  by 
obtaining  exact  solutions  to  an  approximate  equation,  such  as: 

(ao/at)  +  A[Q^,Q^'*'^](aQ/ax)  =  o  (23) 

where,  A[Q^,Q^'‘'^],  is  Roe's  matrix. 

The  constant  matrix,  A[Q^,Q^'*'^],  is  to  be  constructed  so  as 

to  represent  local  conditions.  Roe^®  required  the  matrix  to  have 
the  following  set  of  properties,  called  Property  U  since  it  is 
intended  to  insure  uniform  validity  across  discontinuities: 

( i )  It  constitutes  a  linear  mapping  from  the  vector  space  Q 
to  the  vector  space  F. 

(ii)  As  approaches  Q,  A[Q^,Q^'*‘^]  approaches  A[Q],  where 

A=(aF/aQ). 

(iii)  For  any  Q^,  :  A[Q^,Q^"*’^]*(Q^'''’ ,Q^)  *  F^'*’^  -  F^. 

(iv)  The  eigenvectors  of  A  aiQ  linearly  independent. 

Once  the  matrix  A  is  obtained,  then  its  eigenvalues  are 
the  wavespeeds  of  the  Riemann  prc±)lem  and  the  right  eigenvector, 

r^,  is  proportional  to  the  change  in  d^jendent  variables,  dQ, 
across  each  characteristic  curve  associated  with  a  specific 

eigenvalue  \^.  Th' ref ore  there  are  two  ways  of  obtaining  the 
flux  vector  at  i+1/2. 


■i+1/2 


=  f. 


m 

V 

3=1 


^-(j)  ^(j) 


where  the  (-)  superscript  indicates  that  the  suirmation  includes 
only  those  characteristic  curves  that  h.T.>^  n<^>qativG  slope,  that 
is,  only  those  associated  with  negative  eigenvalues.  Another  way 
is 


(24) 
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(25) 


^i+1/2  =  ^i+1 


m 

V 


j=1  3 


where  the  (+)  superscript  indicates  a  sunration  only  over 
positive  eigenvalues.  An  averaging  of  these  two  flux  vector 
r^aresentations  yields 


^i+1/2  =  +  fi+iJ  -  (26) 

Ekjuation  (26)  may  also  be  written  as 

^i+1/2  =  (V2)[fi  +  fi+i]  -  (1/2)  |A(Q^,Q^^^)1  (Qi+i-Qi)  (27) 

lb  construct  Roe's  matrix,  A[Qj_,Qj^^^] ,  for  the  three- 
dimensicxial  Baler  equations  the  information  on  each  side  of  an 
interface  is  used  to  ocrpute 


p  =  (p^  p^^^ 

(28a) 

(28b) 

(28c) 

(28d) 

H  =  (p^Hi  +  p^^^Hi^,)/(p.  +  p^^^) 

(28e) 

a2  =  (yP/p)  =  (y-l)[H  -  (1/2)(u2  +  v2  +  w2)] 

(28f ) 

where  p  is  the  density,  P  is  the  pressure,  u,  v,  and  w  are  the 
velocity  ccnponents,  and  e  is  the  total  energy 

e  =  (P/y-l  )  +  (p/2)(u2  +  v2  +  w2) 

and  H  is  the  total  enthalpy 

H  =  (l/p)(e4-P). 

Roe  and  Pi)te'2  showed  that  these  averages  cure  the  only  ones  that 
have  the  required  prcperties. 

b.  Higher  Order  Spatial  Accuracy 

The  numerical  flux  derived  thus  far  leads  only  to  a  first 
order  accurate  scheme.  A  family  of  numerical  fluxes  of  high 
order  accurate  TVD  schemes  using  Roe  averaging  was  introduced  by 

Osher  and  ChaJcravarthy.'^  A  slight  variation  of 
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that  numerical  flux  vector  at  cell  face  i+1/2  that  was  found  to 
work  well  in  practice  was  obtained  by  using  all  Roe  variables  and 
metrics  corresponding  to  cell  face  i+1/2  for  the  conputation  of 
the  eigenvalues  and  eigenvectors  given  by 

m  _  (  i  ) 

^i+1/2  ^  ^^^Qi^^i+1/2  '-''fj, i+1/2  ^i+1/2 

m  . 

I^{(1-y)/4  [Lt(-i,i)  -  Lj(3,-1)]  + 
(l+y)/4  [Lt(i,-l)  -  LT(l,3)]>r[^^/] 

where 

±  _  .-(j) 

j,i+P/2  ■  i+1/2  '-■<j,i+p/2 

and 


'-^j,i-1/2 

-  ' i+1/2 

(Oi-Qi.,) 

(30b) 

i+1/2 

-  1  ^ 

'i+1/2 

(30c) 

^^j,i+3/2 

1  (  j  ^ 

'i+1/2 

(30d) 

Lj(i  ,n)  = 

:  niirOTod((x*j 

i+l/2' 

(30e) 

minmod(x,y)  =  sign(x)  max{0,min[| x|  ,y  sign(x)]} 

(30f ) 

b  =  (3- 

» 

>. 

1 

>. 

(30g) 

The  parameter  b  is  a  caipression  parameter.^®  The  variables 

I  ^  ^  \  r^  ^  ^  and  \  ^  ^  ^  are  evaluated  using  the  Roe  averaged 
variables.  The  flux  vectors  f^  =  f(Qj^)  in  the  first  term  on  the 

ric^t  hand  side  of  Equation  (29)  is  evaluated  using  the  dependent 
variable  vector  cis  indicated  (not  the  Roe  averaged  variables), 
but  the  subscript  i+1/2  on  this  bracket  indicates  the  metrics  at 
the  i+1/2  are  to  be  used. 

When  two  arguments  are  used  in  the  minttcd  (MINMOD)  operator, 
the  operator  returns  zero  if  the  arguments  are  of  opposite  sign, 
and  the  smaller  absolute  value  of  the  arguments  otherwise.  The 
consequences  of  this  throughout  the  flow  are  discussed  in 
Reference  14. 

Another  limiter  that  is  used  here  is  the  the  highly 
ccrpressive  flux  limiter  due  to  Roe^*^  called  Superbee  (SUPBEE). 
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+  + 

Lj(i ,n)  =  Lj(n,i  ) 

The  scheme  is  independent  of  y  when  the  Superbee  limiter  is 
enployed.  IMs  would  then  appear  to  be  like  a  second  order  Fromn 
scheme  which  occurs  for  y  =  0. 

c.  Inplicit  Solver 

Tb  siirplify  the  solution  matrix,  linearization  of  only  the 
first  order  numerical  flux  is  considered.  The  argument  list  of 
the  flux  vector  and  Jacobian  matrix  is  e>q>anded  in  this  section 
to  include  the  metrics,  denot«3  by  M. 

A  way  of  formulating  the  flux  difference  can  be  carried  out 
by  using  Equation  (25)  at  cell  face  i+1/2  and  Efjuation  (24)  at  cell 
face  i-l/2.  By  using  a  similar  development  as  presented  by 

Whitfieid^,  Equation  (26)  can  )De  written  as 

The  linearized  flux  difference  at  cell  i  using  Eqs.  (31a)  and 
(31b)  is 


^i+1/2  ^i+l/  ^  "  ^i+1/2  ^i+l/^  ^i+1/2  ^^+1  ^i+1/2  ^^i+l 

+  -  si.,/2  -  +  st.,/2]AQ; 

By  neglecting  the  location  of  the  metrics  in  the  Jacobian 
matrices.  A,  in  Equation  (32),  then  becomes 


^n  +  1  -n  +  1  -yn 

^i+1/2  ^i+1/2  "  ^i+1/2 


f i+1/2 


Equation  ( 33 )  results  in  the  same  fonn  of  the  solution  iratrix  as 
that  ctotained  using  the  flux- vector  split  scheme  discussed 
earlier  (Section  2.2).  The  Jacobian  matrices  in  Equation  (33); 


(31a) 

(31b) 

(32) 

(33) 
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however,  are  the  Roe  matrices,  A  ,  whereas  the  Jacobian  matrices 
in  the  flux-vector  split  scheme  correspond  to  the  true  partials 
of  the  positive  and  negative  flux  vectors.  The  solution  matrix 
resulting  from  Equation  (33),  with  Roe  matrices,  was  used,  as  well  as 
this  sane  solution  matrix  but  with  the  Jacobian  matrices 
corresponding  to  the  flux-vector  split  scheite.  In  all  results 
ctotained  thus  far,  the  Jacobian  matrices  corresponding  to  the 
flux-vector  split  schene  gave  iitproved  convergence  rates  and  a 
more  robust  scheme  than  the  Roe  matrices.  The  behavior  is 

similar  to  what  happens  in  the  flux-vector  splitting,  i.e.  A"*"  =  T/\T“\ 
are  used  in  place  of  the  true  Jacobian  matrices  obtained 
from  differentiating  the  positive  and  negative  flux  vectors. 

Because  the  Jacobian  matrices  corresponding  to  the  flux- 
vector  split  scheme  work  better  on  the  left  hand  side  in  the 
solution  matrix  than  the  Roe  matrices,  this  is  the  scheme 
presently  being  employed  in  the  flow  solver.  This  implementation 
of  the  Roe  scheme  is  particularly  sinple,  in  that  the  computation 

of  the  residual  vector,  R*^,  in  the  flux-vector  split  scheme  need 
only  to  be  r^ laced  with  the  Roe  scheme. 

4.  BOUNDARY  CONDITIONS 

To  derive  the  characteristic  variable  boundary  conditions, 
the  Euler  equations  need  to  be  cast  in  characteristic  variable 
form.  Therefore,  the  characteristic  variables  are  given  first 
frcm  Reference  2,  followed  by  the  derivation  of  the  boundary 
conditions  for  the  specific  cases  of  supersonic  inflow  and 
outflow,  and  impermeable  surfaces. 
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From  Reference  2,  the  characteristic  variables  correspond  in 
order  to  the  eignvalues 


kyU  +  kyV  +  k2W  = 

(34a) 

II 

0k 

(34b) 

II 

0k 

(34c) 

0k  + 

(34d) 

II 

0k  -  c|  Vk| 

(34e) 

where  k=  ^ ,  or  ^  . 

a.  Characteristic  Variable  Boundary  Conditions 

The  boundary  conditions  are  derived  below  assuming  locally 
one-dimensional  flow.  This  assunption  is  prdDably  better  for  far 
field  boundary  conditions  than  for  boundary  conditions  applied  on 
or  near  surfaces.  However,  numerical  es^jeriments  using  zero 
pressure  gradient,  normal  pressure  gradient,  and  locally  one¬ 
dimensional  characteristic  variable  boundary  conditions,  indicate 
that  similar  results  can  be  obtained  using  any  of  these  three 
methods  for  inpermeable  wall  boundary  conditions  as  long  as  the 
grid  is  not  extraordinarily  coarse.  For  the  corputations 
performed  thus  far,  the  locally  one-dimensional  characteristic 
variable  irrpermeable  wall  boundary  conditions  are  to  be  preferred 
over  the  other  methods;  hence,  this  method  is  developed  below. 

The  eigenvalues  indicate  a  direction  in  corrputational  space, 
where  one  particular  eigenvalue  is  associated  with  one  particular 
characteristic  variable.  Each  eigenvalue  indicates  the  direction 
across  the  k  =  constant  computational  surface  that  information 
contained  in  the  associated  characteristic  variable  propagates. 
This  result  is  the  basis  for  determining  the  boundary  conditions 
referred  to  here  as  characteristic  variable  boundary  conditions. 
Boundary  conditions  are  now  developed  for  the  following  five 
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specific  cases: 


1 .  far field-supersonic  inflow 

2.  far  field-supersonic  outflow 

3.  farfield-subsonic  inflow 

4.  farfield-subsonic  outflow 

5.  irrpermsable  surface 

Farfield-Supersonic  Inflow 

This  is  a  situation  where  all  eigenvalues  have  the  same  sign. 
Because  flow  is  coming  into  the  conputational  domain,  all  flow 
variables  are  specified. 

Farfield-Supersonic  Outflow 

This  is  another  situation  where  all  eigenvalues  have  the  same 
sign.  Because  flow  is  leaving  the  conputational  domain  all  flow 
variables  at  the  boundary  must  be  obtained  from  the  solution  in 
the  conputational  domain.  All  flow  variables  are  extrapolated 
from  inside  the  conputational  domain  to  the  boundary. 

Farfield-subsonic  Inflow 

This  situation  is  characterized  by  four  eigenvalues  of  the 
same  sign  and  one  of  differing  sign.  For  the  subsonic  inflow 
case  shown  in  Figure  2a  with  the  flow  in  the  direction  of 
increasing  conputational  coordinate  k,  the  first  four  eigenvalues 
are  positive  and  the  fifth  is  negative.  For  the  subsonic  inflow 
case  shown  in  Figure  2b  with  the  flow  in  the  direction  of 
decreasing  computational  coordinate  k,  the  first  three  and  fifth 
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b.  Flow  In  Dlroctii 


Figure  2.  Ccnputational  Coordinates  Schatia 
Variable  Boundary  Conditions 
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eigenvalues  are  negative  and  the  fourth  eigenvalue  is  positive. 
For  a  totally  general  three-dirnensional  code  either  situation  in 
Figure  2  could  occur.  Solving  for  the  deperdent  variables  across 
the  boundary  yields  (Figure  2): 


=  (l/2){Pa  +  P,  +  sign{xJ.)p^Co[kx(Ua-u,  ) 

+  ky(v^-V|  )  +  kziwg-Wj  )]}  (35a) 

Pb  =  Pa  +  C^^b  -  (35b) 

“b  =  “a  +  kx[(Pa-Pb^/po^olsign(\J^)  (35c) 

=  '^a  + 

^b  =  '^a  + 


Note  that  these  signs  correspond  to  the  sign  of  the  first  three 
eigenvalues,  and  hence  this  is  a  means  of  writing  the  code  for 
general  aEplications  with  arbitrary  orientation  of  the 
oonputation  ooordinates.  The  point  a  is  outside  the 
conputational  domain,  point  b  is  on  the  coirputational  boundary, 
and  point  i  is  inside  the  conputational  domain. 

Farfield-Subsonic  (Outflow 

This  situation  is  also  characterized  by  four  eigenvalues  of 
the  same  sign  and  one  of  opposite  sign.  The  development  of 
subsonic  outflow  boundary  conditions  is  similar  to  that  for 
subsonic  inflow,  and  Figure  (2)  can  be  used  again  for  illustration. 
However,  for  subsonic  outflow  only  one  characteristic  variable  is 
specified  and  four  are  determined  from  information  inside  the 
conputational  domain,  whereas,  for  subsonic  inflow  four 
characteristic  variables  were  specified  and  one  was  determined 
from  information  inside  the  conputational  domain.  Note,  however, 
that  although  the  equations  are  identical,  point  a  is  now  inside 
the  conputational  domain  and  point  b  is  outside  the  conputational 
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dcmain;  whereas,  just  the  opposite  was  true  for  subsonic  inflow. 

The  remaining  four  equations  can  be  solved  for  the  remaining 
four  variables  giving 


Pb  =  Pa 

(36a) 

Pb  =  Pa  ^ 

(36b) 

{36c) 

^b  =  Va  +  ky[(Pa-Pb)/p^c^]sign(\^  ) 

(36d) 

Wb  =  '^a  +  kzC<Pa-Pb^/po^o]sign(\|'  ) 

(36e) 

Inpermeable  Surface 

For  a  boundary  across  which  there  is  no  flow  the  first  three 
eigenvalues  given  by  Equation  (34)  are  zero,  the  fourth  is  positive, 
and  the  fifth  is  negative.  One  condition  must,  therefore,  be 
specified.  The  condition  specified  is  that  there  is  no  flow 
across  the  boundary. 

Finite  volume  codes  only  require  the  pressure  at  an 
iitpermeable  boundary;  however,  to  facilitate  the  handling  of 
points  near  boundaries  and  aid  in  code  vectorization,  frfiantom 
points  are  used  in  the  present  version  of  the  code.  The  use  of 
phantom  points  requires  information  of  variables  other  than 
pressure  to  ensure,  for  example,  zero  flow  across  an  impermeable 
boundary. 

Ihe  solution  for  the  remaining  four  variables  is  es  follows: 


It 

QL4 

Pr 

+  k^Wj.) 

{37a) 

Pr 

+ 

1 — 1 

O' 

1 

1 _ 1 

(37b) 

Ur 

kjj(kj{Ur  +  kyVj.  + 

k2Wr) 

(37c) 

''b  = 

vr 

-  ky(kxUi-  +  kyVj.  + 

^z''^r  ^ 

(37d) 

'^b  = 

'^r 

k^dc^Ur  +  kyVj-  + 

k2Wr) 

(37e) 

where  the  point  r  is  the  center  of  the  fiisf  cell  from  the 
boundary  and  the  minus  sign  in  Equation  (37a)  is  used  if  r  is  in  the 
positive  k  direction  from  the  boundary,  and  the  plus  sign  is  used 
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if  r  is  in  t±e  negative  direction  from  the  boundary. 


b.  Phantom  Points 

Phantom  points  are  denoted  by  the  subscript  p.  The  points 
are  ctotained  from  the  relations 


Pp  =  2Pt, 

-Pin 

(38a) 

Pp  =  ^l^b 

(■'in 

(38b) 

^p  =  2u^j 

-  ^in 

(38c) 

Vp  =  2vj^ 

-  ''in 

(38d) 

''p  =  2Wb 

-''in 

(38e) 

where  the  subscript  "in”  refers  to  the  center  of  the  first  cell 
inside  the  computational  domain  and  can  be  any  of  the  points  a 
or  r  used  in  this  section.  For  example,  the  phantom  points  for 
an  impermeable  surface  are 


=  Pr 

^  2p^Co(kj{Ur  +  kyVr 

+  k2W3.) 

(39a) 

^'p 

=  ^'r 

+  2(Pj,  -  Pj 

J/c§ 

(39b) 

Up 

=  Ur 

2kjj(kj^r 

+  kyVr  + 

(39c) 

''P 

=  Vr 

-  2ky(kj^r 

+  icyVr  + 

(39d) 

Wp 

=  Wr 

-  2k2(kjjUr 

+  kyVr  + 

(39e) 

The  velocity  vector  components  are  the  same  as  those  used  by 
Jaoocks  and  Kneile.^^ 

c .  Blocked  Grids 

Blocked  grids  are  handled  in  much  the  same  manner  in  as  much 
as  phantom  or  image  points  are  required  at  block  boundaries  to 
allow  for  second  order  differencing  across  the  block-to-block 
interfaces.  In  the  present  code,  the  dependent  variables  that 
correspond  to  a  given  block  inter  far:-  ar-^  st>''red  so  that  the 
adjacent  block  will  have  access  to  those  values  when  they  are 
required.  These  interfaces  are  identified  in  the  CUT  (Section 
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4.5)  and  PASTE  (Section  4.24)  subroutines  and  the  values  are 
retrieved  and  stored  at  the  proper  times  ty  the  GETBC  (Section 
4.14)  and  FUTBC  (Section  4.26)  subroutines.  For  this  steady- 
state  code,  the  dependent  variables  used  are  "unsynchronized" 
because  the  code  takes  whatever  time  level  values  are  present 
wh^  they  are  needed.  This  is  ofposed  to  storing  additional  time 
levels  required  for  the  bock-to~block  interfaces  which  would 
cause  an  increase  in  the  amount  of  in-core  storage.  From 
Reference  1 ,  the  difference  between  synchronized  and 
unsynchronized  bDundary  conditions  is  substantial  for  unsteady 
(or  time  accurate)  solutions;  however,  for  the  steady  Ccise  where 
the  error  approaches  zero  there  appears  to  be  no  degradation  in 
convergence  or  stability. 

d.  Zero  Pressure  Gradient  (ZPG)  Ctondition 

For  several  instcinces  a  "starting"  solution  may  need  to  be 
Implemented  for  the  flow  solver  using  a  zero  pressure  gradient 
boundary  condition.  This  type  of  BC  is  used  to  enable  the  code 
to  transition  from  the  freestream  conditions  at  step  "0"  to  the 
characteristic  variable  BCs  (with  inpentieable  wall  BC)  at  step 
NZPGf  1 .  Ihis  allows  the  code  to  "build  up"  to  a  complete  solid 
body  boundary  in  a  specified  nunber  of  zero  pressure  gradient 
iteraticais  (NZPG).  This  BC  forces  the  code  to  allow  some  flow 
through  the  body  for  several  (N2TC)  iterations  as  shown  in  Figure 
(3).  lb  further  enhance  the  codes  ability  to  smoothly  transition 
through  the  beginning  iterations,  the  inputs  allow  the  user  to 
gradually  apply  the  zero  pressure  gradient  boundary  condition 
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along  iirpermeable  surfaces  (NGRAD).  These  boundary  conditions 
are  irtplemented  in  the  following  manner  on  a  constant  k 
infjentieable  surface: 


= 

Pin 

(40a) 

Pp  = 

Pin 

(40b) 

Up  = 

‘^in  -  Pz 

(40c) 

Vp  = 

''in  -  Pz  Vin  ky 

(40d) 

Wp  = 

''in  "  Pz  ''in  ^z 

(40e) 

1"  2.0 

for  zero  pressure  gradient  BCs 

1  2.0  * 

(st^  #/N®AD)  [less  than  or  equal  to  2.0] 

L 

for  gradually  applied  2PG  BCs 
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SECTION  III 


FLOW  SOLVER  INPUT 

Namelist  input  allows  increased  flexibility  in  the  ability 
to  group  variables  and  to  read  those  variables  from  an  input 
deck.  This  capability  has  thus  been  iJtpletnented  in  this  code  to 
allow  a  better  understanding  of  the  interaction  among  the  input 
variables.  The  following  constitute  the  namelist  headings: 

Namelist  /FINPOT/  :  general  flow  solver  input 
Namelist  /ROTATE/  :  transformation  angles 
Namelist  /SURFACE/  :  surface  identification 
Namelist  /BCIN/  :  boundary  corriitions 

1 .  GENERAL  INPOT  /ETNPUT/ 

As  stated  previously,  namelist  FINFUT  makes  up  the  general 
flow  solver  input.  FINPOT  is  the  first  of  the  namelists  to  be 
reed  and  is  only  read  once  at  the  beginning  of  the  main  routine 
(MISSE). 


a.  CFL  (real):  Courant  Nurrber 
The  CFL  variable  is  input  to  set  the  maximum  time  step 
allowed  during  local  time  stepping.  The  default  value  is  <15.0> 
but  is  usually  only  set  that  high  for  siirple  cases  such  as 
airfoils  at  subsonic  conditions.  Typically,  CFL  should  fall 
between  5.0  and  10.0  d^jending  on  angle  of  attack,  bluntness 
ratio  of  configuration  and  Mach  nunber. 

Example:  E$  FINFUT  CFL=8.0,... 
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b.  FSMACH  (real):  freestream  Mach  number 

FSM^CH  sets  the  ratio  of  freestream  velocity  to  local  speed 
of  sound.  Ihe  FSMACH  is  used  to  set  the  initial  conditions  for 
the  velocity  coipanents  in  the  x,  y  and  z  directions  (u,  v,  and 
w).  The  default  value  for  FSMACH  is  <0.95>  but  can  be  set  to 
Mach  numbers  in  the  cotipressible  subsonic,  transonic  and 
supersonic  range. 

Example:  E$  FINPUT  CFL=  8.0,  FSMACH  =  0.80, — 

c.  LIFIAX  (integer):  Lift  axis 

The  LIFTAX  variable  is  input  to  tell  the  code  in  which 
direction  to  place  the  lift  vector.  This  value  is  set  before  any 
rotation  is  made  to  the  right-harbted  coordmate  system  (namelist 
ROTATE);  therefore,  the  LIFIAX  variable  is  set  to  the  positive  y 
direction  (+2),  positive  z  direction  (+3),  the  negative  y 
direction  (-2)  or  to  the  negative  z  direction  (-3).  The  default 
value  is  <+2>  and  a  representation  of  these  situations  is  shown 
in  Figure  4. 

Example:  E$  FINRJr  CFL=8.0,  FSMACH=0.80,  LIFTAX=-2... 

d.  JFRE2  (integer);  frequency  of  updates  to  t)ie  flux  Jacdbians 

The  JFMQ  variable  sets  the  number  of  iterations  at  which 

the  flux  Jacobian  matrices  are  updated.  If  JFRBQ  is  set  to  10, 
the  flux  Jacobian  matrices  are  rejalculate.l  ‘-'■.•ery  10  iterations. 

The  default  valve  is  < 1 >  which  requires  the  matrices  to  be 
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updated  every  iteration  during  the  entire  job.  For  steady  state 
solutions,  as  we  have  here,  results  have  shown  that  the  two 
extremes  {JFBEQ=0  and  JFRE2=1 )  yield  the  best  results. 

Typically/  for  those  cases  with  sharp  gradients  in  the  flow  field 
the  latter  (JFRBQ=1 )  situation  will  enable  the  solution  to 
converge,  although  the  runtime  will  be  increased  dramatically 
(updating  the  flux  Jaccijian  matrices  requires  almost  50%  of  each 
iterations  CHJ  time ) .  For  sinrpler  conditions ,  the  former  input 
(jniE3=0)  will  suffice  and  greatly  reduce  CPU  time. 

Example:  E$  FINPOT  CFL=8.0,  FSMACH=0.80,  LIPTAX=-2, 

jniEQ=5.0, . . . 

e.  reiNT  (integer):  Print-out  specification 

The  PRIwr  variable  sets  the  type  of  print-out  the  user 
requires.  Tne  following  choices  are  made  available  to  the  user: 

PRIwr=0;  forces  and  moments  only 

PRINr=l;  X  versus  Cp  and  forces  and  moments 

PRINr=l23;  x,y,z  versus  Cp  and  forces  and  moments 

The  default  value  for  PRINT  is  <0>. 

Example:  E$  FINPUT  CFL=8.0,  FSMACH=080,  LIFrAX=-2,  JFRBa=50, 

PRINT=  1 , . . . 

f .  NPR  ( integer ) :  number  of  printouts  ; 

NIT  (integer):  number  of  iterations 

The  NPR  variable  sets  the  number  of  prints  the  user 
requires.  NPR,  in  conjunction  with  NIT,  also  sets  the  total 
number  of  iterations  (NPR*NIT).  The  NIT  variable  sets  the  number 
of  iterations  at  which  a  printout  will  occur.  Therefore,  if  the 
user  requires  4  prints  and  a  total  of  200  iterations,  NPR  would 
be  set  to  4  and  NIT  to  50  (NPR=4,  NIT=50).  In  this  case,  every 
50  iterations  a  printout  will  be  accomplished  under  the  format 
set  by  PRINT  ( Section  3.1)  There  are  no  default  values  for 
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either  NPR  or  NTT. 


Exaitple:  E$  FINPUT  CFL=8.0,  F»1ACH=0.80,  LIFrAX=-2,  JFREQ=50 

PRINT=1,  NPR=6,  NIT=100,...$ 

g-  NCYC  (integer):  nuirber  of  cycles  for  residual  printout. 
The  NCYC  variable  sets  the  frequency  of  residual  prints  to 
allow  the  user  to  judge  the  convergence  level  of  his/her  run. 
Residual  printouts  are  formatted  as  follows: 


CYCLE  BLOCK  RIMAX  I  J  K  RTRMS  BIMAX  FIRMS  L2N0RM 


where 


CYCLE:  iteration  number 

BLOCK:  block  number 

Rl^ftX:  maximum  residual  for  density  (R) 

I,J,K:  I,J,K  location  at  which  RTMAX  occurs 
RTRMS:  RMS  value  of  density  (log) 

ETMAX:  Maximum  residual  for  energy  (E) 

FIRMS:  RMS  value  of  energy  (log) 

L2N0BM:  L2  norm  value  for  d^)endent  veiriables  (log) 
NSUP:  Nuntoer  of  supersonic  points  in  entire  field 


An  NCYC  value  of  5  will  yield  a  residual  print  out  every  fifth 
iteration  (cycle).  The  default  value  is  <1>. 


Example:  E$  FINHIT  CFL=8.0,FSMACH=0.80,  LIFrAX=-2, 

JFRBQ=50,  FRINT=1,  NPR=6,  NIT=  100, 
NCYC=10,... 


h.  NSURF  (integer):  number  of  surfaces 
The  NSURF  variable  tells  the  code  how  many  surfaces  are  used 
in  calculating  the  forces  and  moments  for  the  configuration.  The 
default  value  is  < 1 >  and  more  details  will  be  discussed  in 
Section  3.3  (Namelist  /SURFACE/). 

Example:  E$  FINPUT  CFL=  8.0,  FSMACH=  0.80,  LIFTAX=-2, 

JFE^=50,  PRINT=1,  NPR=6,  NIT=100, 
NCYC=10,  NSURF=3,... 


NSUP 
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i.  SPLIT  (integer):  type  of  splitting  technique  enployed. 

Hie  SPLIT  variable  sets  the  type  of  splitting  technique  used 

in  the  code.  Hie  user  has  the  following  choices: 

SPLIT=1;  flux- vector-splitting  of  Steger  and  WarTning 
SPLIT=2;  flux-difference-splitting  of  Roe 

Ihe  default  value  is  Roe's  flux-difference-splitting  <2>  which 
is  less  dissipative  than  flux-vec±or-splitting.  When  SPLIT=2  is 
enployed, a  second  variable,  LIMIT  (Section  3.1 ),  must  be  invoked. 

The  Steger  and  Warming  splitting  technique  (SPLIT=  1 )  is 
available  due  to  its  tendency  to  reach  a  converged  solution  at  a 
faster  rate,  when  this  splitting  technique  is  enployed  the  LIMIT 
input  is  ignored. 

Exanple:  E$  FINPOT  CFL=8.0,  FSMACH=0.80,  LIFI!AX=-2,  JFRBQ=50, 

Hmfr=1,  NFR=6,  NIT=100,  NCYC=10,  NSURF=3, 
SPLIT=2, . . . 

j .  LIMIT  ( integer ) :  type  of  flux  limiter  used  during  Roe  averaging 
The  LIMIT  veuriable  is  used  to  specify  the  type  of  flux 

limiter  used  during  a  Roe  averaged  flux-difference-split  run 
(SPLIT=2).  The  following  choices  are  given  to  the  user: 

LIMIT=1;  calls  subroutine  MINMDD  (Sec  4.22) 

LIMIT=2;  calls  subroutine  SUPBEE  (Sec  4.39) 

For  more  information  on  the  limiters  enployed,  refer  to  the 
sections  describing  the  specific  subroutines.  The  default  value 
is  <1>.  If  SPLIT=1  is  enployed  then  the  LIMIT  input  is  sinply 
ignored  and  no  error  message  is  sent  by  the  code. 

Exanple:  E$  FINPOT  CFL=8.0,  FSMACH=0.80,  LIFTAX=-2,  JFRB3=50, 

PRINT=1,  NPR=6,  NIT=100,  NC:YC=10,  NSURF=3, 
SPLIT=2,  LIMIT=2,  ... 
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k.  BORDER  (integer):  nunter  of  phantcm  points  used 

The  BORDER  variable  is  used  to  set  the  number  of  piiantom 
points  (Section  2.3)  used  in  the  block  boundary  calculations. 

The  following  choices  are  given  to  the  user: 

BORE®R=  1;  first  order  block  boundary  (one  phantcm  point) 

BORDER=  2;  second  order  block  boundary  (two  phantom  points) 

The  default  value  is  <2>  and  there  has  been  no  reason  to  go  first 
order  in  any  test  cases  to  the  present  time. 

Exanple:  E$  FINPOT  CFL=8.0,  PSMACH=0.80,  LIFEAX=-2,  JFREQ=50 

PRINr=1,  NPR=6,  NIT=100,  NCYC=10, 

NSURF=3,  SPLIT=2,  LIMIT=2,  B0RDER=2,... 

l.  GAM  (real):  ratio  of  specific  heats 

The  GAM  variable  allows  the  user  to  inpxit  a  ganma  value  for 
any  type  of  fluid  medium.  The  default  value,  <1.4>,  is  set  for 
air. 

Exanple:  E$  FimUT  CFL=8.0,  FSMACH=0.80,  LIPrAX=-2,  JFRBQ=50 

PRDfr=1,  NFR=6,  NIT=100,  NCyC=10, 

NSURF=3,  SPLIT=2,  LIMIT=2,  B0RDER=2, 
GAM=1.31,... 

m.  RESTRT  (integer):  restart  flag 

The  RESTRT  variable  tells  the  code  whether  or  not  a  restart 
file  is  being  read  in  and  from  what  file  tape  number  the  restart 
file  is  being  read  from.  The  following  delineates  the  users’ 
cptions: 

RESTRT  =  0;  no  restart 

RESTRT  >  1 ;  restart  will  be  read  in  on  unit  RESTRT  and  a  new 
file  will  be  written  out  on  RESTRT  +  1 . 

The  default  value  of  RESTRT  is  <0>. 

Exanple:  E$  FINPOT  CFL=8.0,  FSMACH=0.80,  LIFrAX=-2,  JFREQ=50 

PRINT=1,  NPR=r3,  NIT-100,  NHi’C-lO, 

NSURF=3,  SPLIT=2,  LIMIT=2, 

B0RDER=2,  GAM=1.131,  RESTRT=8,... 


n.  NZPG  (integer):  number  of  zero  pressure  gradient  BC  iterations 
Ihe  NZPG  variable  is  input  to  tell  the  cxxJe  how  many  zero 

pressure  gradient  Loundary  (Condition  iterations  are  to  be  allied 
in  a  specific  job.  The  default  value  is  <10>;  however,  this 
value  is  typically  only  good  for  airfoils  at  fairly  simple  flow 
conditicais  (subsonic  at  low  angles  of  attack).  NZPG  values  of  50 
to  100  are  generally  applied  to  three  dimensional  configurations 
depending  on  Mach  number,  bluntn^s  ratio  of  shape  and  angle  of 
attack. 

Exanple:  E$  FINPOT  CFL=a.O,  Fa^CH=0.80,  LIPEAX=-2,  JFRBQ=50, 
PRINT=1,  NPR=6,  N1T=100,  NCYC=10,  NSURF=3, 
SPLIT=2,  LIMIT=2,  B0RDER=2,  G?M=1.131, 

RESTRT=8,  NZPG=75,... 

o.  MGRAD  (integer):  number  of  gradually  applied  ZPG  BC  iterations 
The  NGRAD  variable  sets  tte  nuntoer  of  iterations  that  zero 

pressure  gradient  boundary  conditions  will  be  applied  gradually 
(Section  2.4).  NGRAD  must  always  be  less  than  NZPG.  The  default 
for  NGRAD  is  <0>  although  typically,  for  tJuree-dimensional 
configurations,  NGE1AD  should  be  from  10  to  50  or  more  depending 
on  angle  of  attack,  Mach  number  and  bluntness  ratio. 

Exarrple:  E$  FINPUT  CFL=8.0,  FSMACH=0.80,  LIFIAX=-2,  Jn<BQ=50, 
PRINT=1,  NPR=6,  NIT=100,  NC:YC=10,  NSURF=3, 
SPLIT=2,  LIMIT=2,  B0RDER=2,  GAM=1.131, 

RESTRT=8,  NGRAD=25,  NZPG=75,... 

p.  RORDEK  (integer):  Order  of  Roe  scheme  approximation 
The  RORDER  input  specifies  the  order  of  approximation  used 

in  the  Roe  averaging.  The  following  chKoices  are  given  to  the 
user: 
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RiCM)ER  =  1  ;  first-order  approxiitation  (no  limiters) 

RCRDER  =  2  ;  second-order  approximation 

RORDER  =  3  ;  third-order  approximation  (MINMOD  only) 

The  default  value  is  <2>  for  a  second-order  approximation. 

Exanple:  E$  FJNHJT  CFL=8.0,  FSMAa)=0.80,  LIFI!AX=-2,  JFRBQ=50, 

PRINT=1,  NPR=6,  NIT=100,  NC:YC=10,  NSURF=3, 
SPLIT=2,  LIMIT=2,  B0RDER=2,  GAM=1.131, 
RESTRT=8,  NGRAD=25,  NZPG=75,  R0RDER=1,  ... 

2.  TRANSPORMATIW  INPUT  (ROTATE) 

The  transformation  input  sec±ion  allows  the  user  to 
manipulate  a  set  of  transformation  matrices  to  rotate  the 
configuration  of  interest  about  any  of  its  three  axes.  This  will 
effectively  give  the  configuration  an  angle  of  attack,  (;y ,  roll 
angle,  ,  or  side-slip  angle,  0 ,  or  any  combination  of 
these  angles  depending  on  how  the  configuration  is  oriented.  The 
orientation  is  input  as  LIPTAX  in  the  set  of  general  input 
(Section  3.1 )  and  tells  the  code  in  which  direction  to  place  the 
lift  vector. 

a.  TRANS  (integers):  order  in  which  the  rotations  will  occur 
TRANS  acts  to  tell  the  code  in  which  order  to  rotate  the 
"Euler"  angles.  The  integer  "1"  will  rotate  the  configuration 
about  the  "x"  axis,  "2"  will  rotate  the  configuration  aJxut  the 
"y"  axis,  and  "3"  will  rotate  the  configuration  alx>ut  the  "z" 
axis.  All  angles  are  counter-clockwise  positive  (Figure  5) 
looking  towards  the  origin.  The  default  values  are  all  set  to 
zero,  effectively  yielding  no  rotations.  If  only  one  rotation  is 
required,  say  about  the  "y"  axis,  then  only  a  "2"  will  be 
entered. 

Exanple:  E$  ROTATE  TRANS  =  2,  PHIY  =  -4.00  $ 
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In  this  case,  the  configuration  will  be  only  rotated  clocJcwise 
about  its  "y"  axis.  The  other  two  angles  are  assumed  to  remain 
zero  and  so  their  values  for  'IKANS  need  not  be  entered.  The  same 
situation  holds  when  only  two  rotations  are  requested. 

Exaitple:  E$  ROTATE  TRANS  =  3,1,  PHIX  =  -2.43  ,  PHIZ  =  5.84$ 

In  this  case,  the  configuartion  will  be  rotated  counter-clockwise 
first  about  the  "z"  axis  and  then  rotated  clockwise  about  the  "x" 
axis.  Note  that  the  order  of  the  actual  angles  (PHIX  and  PHIZ) 
is  inconsequaitial. 

It  is  possible  to  rotate  abcaat  the  same  axis  more  than  once 
in  a  given  run;  however,  this  case  is  limited  by  the  fact  that 
the  magnitude  of  rotation  must  remain  the  same  for  both  rotations. 

Example:  E$  ROTATE  TRANS  =  3,1,3,  PHIX=  2.43  ,  PHIZ=  5.84 

Here  the  user  is  first  rotating  about  the  ”z"  axis  in  a  counter¬ 
clockwise  direction,  then  about  the  "x"  axis  in  a  counter¬ 
clockwise  direction  and  then,  again,  about  the  "z"  axis.  Both 
rotations  about  the  "z"  axis  are  of  the  same  magnitude  and 
direction.  There  is  no  logic  in  the  code  to  allow  the  user  to 
rotate  about  the  same  axis  more  than  once  with  different  "PHI" 
values. 

b.  PHIX  (real):  rotation  about  the  "x"  axis  (degrees) 

PHIY  (real):  rotation  about  (In--  "y"  axis  (tiegrees) 

PHIZ  (real):  rotation  about  the  "z"  axis  (degrees) 
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The  "PHI"  angles  are  sirtply  the  amount  of  rotation  about  the  given 
axis.  Again,  the  angles  are  given  in  the  counter-clockwise 
positive  direction  and  work  in  conjunction  with  the  IPANS  input 
(Section  3.2)  and  the  LIFIAX  input  (Section  3.1).  Ihese  angles 
must  be  input  in  degrees  as  stated  above  and  are  converted  into 
radians  within  the  code. 

Exanple:  E$  ROTATE  TOANS=  1,  PHIZ=  -12.0 

In  this  exanple,  the  user  is  attenpting  to  attain  an  angle  of 
attack  of  12.0  degrees.  The  LIFIAX  value  entered  in  the  general 
input  section  (FINRJT)  was  specified  as  <2>,  the  positive  "y" 
axis.  This  means  that  to  obtain  a  tail  down,  positive  angle  of 
attack  of  12  degrees,  the  configuration  is  rotated  clocdcwise 
about  the  "z”  axis  effectively  yielding  the  required  angle. 

3.  SURFACE  SPECIFICATION  (SURFACE)) 

Ihis  set  of  inputs  allows  the  user  to  delineate  the  surfaces 
he/she  wants  to  calculate  forces  and  mcments  on  or  about  and  sets 
the  reference  lengths  in  the  "x",  "y”  and  "z"  directions  as  well 
as  the  total  reference  length  and  area. 

a.  Surface  Name  (sixteen  alpha-numeric  characters) 

A  surface  name  (maximum  of  sixteen  characters)  is  input 
prior  to  every  namelist  SURFACE.  This  aids  in  identifying  the 
surfaces  that  will  be  sp)ecified  in  the  force  and  moment  print-out 
section. 
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b.  SREF  (integer):  surface  number 

The  value  for  SREF  is  siitply  a  numerical  identifier  for  the 
code  and  corresponds  to  the  SURF  values  to  be  input  with  the 
namelist  BCIN  cards  in  Section  3.4.  SREF=  1  must  always  be  the 
entire  configuration  required  in  the  force  and  mcment 
calculations. 

c.  XREF  (real):  reference  displacement  off  the  "x"  axis 
YREF  (real):  reference  displacement  off  the  "y"  axis 
ZflEF  (real):  reference  displacement  off  the  "2"  axis 

The  value  for  XREF,  YREF  and  ZREF  allows  the  user  to  set  a 
position  about  which  the  code  will  calculate  moments.  Default 
values  are  listed  below: 

XREF  =  0.0,  YREF  =  0.0,  ZREF  =  0.0, 

This  shows  that  the  default  values  for  XREF,  YREF  and  Z(?EF  lie  on 
the  "x",  "y”  and  "z”  axes  respectively. 

d.  AREF  (real):  reference  area 
LREF  (real):  reference  length 

The  AREF  and  LREF  values  are  input  so  as  to  allow  the  user  to 
set  a  reference  area  (planform,  cross-sectional  or  any  sectional 
area  required)  and  a  reference  length,  be  it  chord  length,  axial 
bocty  length  or  any  section  length  required.  The  default  values 
for  these  input  values  are  as  follows: 

AREF  =  1.0,  LREF  =  1.0. 

An  exanple  of  the  Surface  specification  input  is  given  as  follows 


Exanple:  Generic  body 

E$  SURFACE  SREF=  1,  AREF=  17.93  ,  LREF=  142.53,  XREF=  3.2 
Fin  leading  edge 

E$  SURFACE  SREF=  2,  AREF=  56.78  ,  LREF=  3.42,  XREF=  .855 


Nose  cone 

E$  SURFACE  SREF=  20,  AREF=  23.45 

Ibis  exanple  shows  the  surface  name  of  the  first  section  of  body 
to  be  specified,  "Generic  body",  and  it  must  enconpass  the  entire 
configuration.  Ihis  is  labeled  as  SREF=  1  (corresponding  to  all 
of  the  sections  specified  in  the  BCIN  input  with  SURF=  1 )  with  a 
cross-sectional  area  of  17.93  units,  a  reference  length  of  142.53 
units  and  a  +3.2  unit  displacement  on  the  "x"  axis  for  the  mcment 
reference  point.  The  second  section  (corresponding  to  all 
sections  of  BCTYPE*  1  that  includes  a  2  in  the  SURF  card)  to  be 
specified  is  the  "Fin  leading  edge"  with  a  fin  planform  area  set 
at  56.78  units,  a  chord  length  of  3.42  units  and  the  moment 
reference  point  set  at  the  fin  quarter  chord.  There  can  be  a 
maximum  of  20  surfaces  specified,  as  indicated.  All  input  values 
cire  echoed  out  at  the  end  of  a  printout  along  with  the  forces 
and  momaits  calculated  for  each  specified  surface. 

4.  BOUNDARY  CONDITIONS  (BCIN) 

The  following  set  of  inputs  requires  the  user  to  specify  a 
boundcury  condition  type  for  all  points  on  the  boundary  of  a  given 
ccrputational  grid.  There  will  be  a  total  of  six  ccnplete  sides 
in  a  three-dimensional  grid  and  four  sides  in  a  two-dimensional 
grid. 
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a.  BCTYPE  (integer):  boundary  condition  type 

Itie  following  shows  the  boundary  condition  types  to  specify 
for  each  patch  of  grid: 

BCTYPE  =  0;  terminates  input 
This  BCIYPE  is  the  flag  that  tells  the  flow  solver  that 
there  are  no  more  boundary  condition  cards  to  be  read.  This  card 
should  be  the  last  line  of  the  input  stream. 

BCTYPE  =  1;  itrpermeable  surface  (solid  body) 

This  BCTYPE  sets  the  solid  body  boundary  patches  for  the 
code  and  flags  the  main  routine  to  call  the  impermeable  surface 
boundary  condition  subroutine  (SDLID-Section  4.37).  The  SURF 
input  (Section  3.4)  is  required  for  BCTYPE  =  1. 

BCTYPE  =2;  far field  conditions 

This  BCTYPE  defines  the  patches  of  grid  on  which  the  code 
iitplements  the  inflow/outflow  boundary  conditions  discussed 
earlier  (Section  2.4)  and  flags  the  code  to  call  the  farfield 
boundary  condition  subroutine  (FftRFLD-Section  4.10) 

BCTYPE  =  3;  synmetry  case  (reflection  plane) 

This  BCTYPE  allows  the  user  to  define  a  plane  of  synmetry 
for  use  in  axisynmetric  cases  (Figure  6). 

BCTYPE  =  4;  block-to-block  information 

This  BCTYPE  defines  the  block-to-blodc  interfaces  on  vhich 
the  flow  solver  must  specify  image  (phantom)  points  (Figure  7). 
For  this  BCTYPE,  both  BLKA  and  BLKB  inputs  are  required. 

BCTYPE  =  5;  block-to-block  across  a  singularity 

This  BCTYPE  specifies  the  boundary  condition  across  a 
singularity,  such  as  the  stagnation  line  out  the  nose  of  an 
aircraft  (Figure  8),  to  give  an  image  point  position  for  use  in 
differencing  across  the  zero  area  cells  created  by  the 
singularity.  For  this  BCTYPE,  both  BLKA  and  BLKB  inputs  are 
required. 

BCTYPE  =  6;  singularity  on  a  reflection  plane 
This  BCTYPE  specifies  the  boundcury  condition  across  a 
singularity,  such  as  the  stagnation  line  out  the  nose  of  an 
aircraft,  when  that  singularity  lies  on  a  reflection  plane 
(Figure  9).  This  input  allows  the  code  to  go  first  order  in  the 
first  cell  along  the  stagnation  line  since  there  is  no  cell  to  be 
specified  as  the  image  position  for  image  point  information. 

BCTYPE  =7;  blockout  conditions  (no  calculations) 

This  BCTYPE  allows  the  user  to  specify  a  portion  of  grid  in 
which  to  make  no  flow  solver  calculations .  The  code  simply  skips 
over  this  region  as  if  it  did  not  exist;  however,  the  rest  of  the 
BCIN  cards  must  reflect  this  change. 

Any  boundary  condition  type  can  be  specified  on  any  patch  of 
grid  by  stating  the  type  (BCTYPE)  of  boundary  condition  and  the 
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SIDE  VIEW 


Figure  9.  Singularity  on  a  Reflection  Plane  (BCTYPE  =  6) 
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start  (STARTA)  and  end  (ENDA)  points  for  a  constant  plane  in  any 
block  (ELKA)  of  corputational  grid. 

b.  SURF  (integers):  surface  numbers  with  BCTYPE  =  1 
This  SURF  input  allows  the  user  to  indicate  in  which 

surfaces  this  patch  of  impermeable  surface  (BCTiTE  =  1 )  is 
included.  As  indicated,  this  input  is  only  valid  for  a  BCIYPE 
equal  to  1  since  only  solid  body  boundaries  are  of  interest  for 
force  and  moment  calculations. 

Exanple:  E$  BCIN  BCrYPE=  1  ,  SURF=  1,  2,  20,  BLKA=  ...  , 
STARTA= . . . ,  ENDA=  . . . 

In  this  example,  surfaces  numbered  1,  2  and  20  fron  the  namelist 
SURFACE  input  (SREF=1,  SREF=2  and  SREF=20)  all  have  this 
particular  patch  of  solid  body  grid  included  in  their  force  and 
moment  calculations. 

c.  BLKA  (integer  ):  block  number  for  a  patch  of  grid 
STARTA  (integers):  start  location  for  a  patch  of  grid 
ENDA  (integers):  end  location  for  a  patch  of  grid 

This  set  of  inputs  specifies  the  location  of  a  patch  of 
corputational  grid  to  be  classified  as  a  particulcir  type  of 
boundary  condition  with  the  BCTYPE  iiput.  BLKA  lets  the  code 
know  which  block  the  patch  is  in,  STARTA  lets  the  code  know  the 
starting  index  of  the  patch  and  ENDA  specifies  the  ending 
location  index  for  the  patch.  The  patch  must  specify  a  constant 
corputational  plane,  and  therefore,  one  irdex  must  remain 
constant  in  the  STARTA  and  ENDA  input  corbination. 

Exanple:  E$  BCIN  BCrYPE=  2,  BLKA=  2,  STARTA=  1,1,1, 

ENDA=  100,25,1 
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This  esxaitple  is  specifying  a  constant  K  plane  as  a  farfield 
bDundary  in  block  2.  Note  that  both  the  I  and  J  indices  are  in 
increasing  order  from  start (STARTA)  to  end(ENDA)  and  that  the  K 
index  is  held  constant. 

d.  BLKB  (integer  ):  block  number  in  block-to-block  conditions 
STARTB  (integers):  start  location  for  block-to-block  conditions 
ENEB  (integers):  end  location  for  block-to-block  conditions 
The  ”B"  inputs  are  only  required  when  a  block-to-block 
boundary  condition  is  ntplemented  (BCIYPE  equals  4  or  5).  This 
set  of  inputs  tells  the  code  which  patch  of  grid  in  the  first 
block  (STARTA,  ENDA,  BLKA)  corresponds  to  which  patch  of  grid  in 
the  second  block  (STARTB,  ENDB,  BLKB).  Using  this  input,  the 
code  passes  boundary  information  between  blocks  A  and  B. 

Ebcanple:  E$  BCIN  BCTYPE=  4, 

BLKA=  3,  STARTA=  1  ,1  ,1,  ENDA=  100,  25,  1, 

BLKB=  8,  STARrB=l25,75,1,  EMB=  25,  50,  1 

In  this  example,  the  user  has  specified  that  block-to-block 
conditions  exist  between  blocks  3  (BLKA)  and  8  (BLKB).  This 
occurs  at  a  constant  K  plane  (as  irdicated  by  the  third  index  in 
both  the  "A"  and  "B"  inputs  being  held  constant).  Note  that  the 
I  and  J  indices  in  the  "A”  input  are  in  increasing  order  fran 
start  to  end,  but  that  this  is  not  required  for  the  "B"  input. 
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SBcri^i  IV 


FLCW  SOLVER  DESCRIPTION 


The  Program  EAGLE  -  Flow  Solver  consists  of  a  main  program 
(MISSE),  several  subroutines  (39)  and  accesses  system  routines 
resident  on  the  Cray  X-MP  (X-MP  version  only).  The  Flow  Solver  is 
specifically  designed  to  execute  on  a  carpiter  with  a  high  speed 
secondary  memory  such  as  the  solid  state  dislc  (SSD)  device  on  a 
Cray  X-MP.  However,  a  Cray  2  version  of  this  code  exists  that 
makes  this  requirement  obsolete.  An  extensive  software 
engineering  effort  has  been  made  to  enhance  the  codes  ability  to 
efficiently  store  and  transfer  data  for  use  within  the  flow 
solver.  Ribbon  vector  storage  with  dynamic  memory  management  is 
enployed  to  prevent  wasting  carputer  memory  when  working  with 
blocks  of  widely  varying  dimensions.  Only  one  block  is  in  memory 
at  a  tine  while  all  other  blodcs  are  stored  on  SSD  [X-MP  version 
only].  The  memory  required  to  field  all  arrays  associated  with  a 
particular  block  is  determined  by  the  code  and  the  memory  size 
will  be  increased  or  decreased  as  each  block  is  operated  on  in 
turn. 

The  term  software  engineering  is  used  loosely  fiere  to 
eitphasize  the  fact  that  considerable  concern  has  been  given  to 
the  efficient  use  of  the  fiardware  be  utilized.  "Four 
properties  that  are  sufficiently  general  to  be  accepted  as  goals 
for  the  entire  discipline  of  software  engineering  are 
modifiability,  efficiency,  reliability  and  understandability"^®^ 
These  essential  tenets  hiave  been  kept  in  mind 
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throuc^xjut  the  develcpnent  of  the  Program  EAGLE  -  Flow  Solver. 

1 .  MAIN  ROUTINE  (MISSE) 

The  main  program  in  the  flow  solver  is  primarily  responsible 
for  the  initialization  of  the  d^sendent  variables  (Section  2.1 ), 
initialization  of  the  ribbon  vector  storage,  the  dynamic  memory 
allocation,  the  reading  of  namelist  inputs  (Section  3. 1-3.4),  the 
reading  and  writing  of  data  files  (ccrputational  grid  or  restart 
files),  and  the  prefer  sequencing  of  subroutine  calls.  In 
addition,  the  MISSE  routine  is  required  to  set  defaults  for 
namelist  iiput  values,  check  values  that  are  set  ty  namelist  for 
any  possible  incorrect  input  and  to  ensure  that  enough  menory  is 
allocated  at  any  given  time  during  the  running  of  the  code. 

Before  the  main  iteration  loop  the  following  subroutines  are 
called: 

(la)  CALL  MEMORY  =>  gives  MEMORY  the  username  from  job 
control  list  (JCL).  [X-MP  version  only] 

(lb)  CALL  MEMORY  =>  gives  MEMORY  users  requested  amount  of 
memory  from  JCL  and  gets  back  the  maximum  allowable  memory 
available  on  the  system.  [X-MP  version  only] 

(lc)  CALL  MEMORY  =>  gives  MEMORY  the  current  amount  of 
memory  being  requested  Iqy  the  Flow  Solver  and  gets  back  the 
minimum  initial  amount  available.  [X-MP  version  only] 

(2a)  CALL  RCHX  =>  gives  ROTX  the  magnitude  of  the  requested 
rotation  after  having  read  the  ROTATE  input.  MISSE  gets  back  an 
updated  transformation  matrix  through  OOW1GN  /TRANS(3,3)/. 

(2b)  CALL  ROTY  =>  gives  ROTY  the  magnitude  of  the  requested 
rotation  after  having  read  the  ROTATE  input.  MISSE  gets  loack  an 
updated  transformation  matrix  through  OOM4CN  /'IRANS(3,3)/. 

(2c)  CALL  ROTZ  =>  gives  ROIZ  the  magnitude  of  the  requested 
rotation  after  having  read  the  ROTATE  input.  MISSE  gets  back  an 
updated  transformation  matrix  through  COMMON  /TRANS(3,3)/. 

(3)  CALL  SETMEM  =>  gives  the  current  location  of  the  last 
memory  location  in  the  ribbon  vector  and  allocates  the  proper 
amount  of  memory.  [X-MP  version  only] 
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(4)  CALL  CHECK  =>  gives  CHEXZK  all  bouncJary  cxjndition  input 
information  after  having  read  a  BCIN  input  card.  MISSE  gets  back 
a  thorough  check  of  all  boundary  condition  input  ensuring  that 
all  locations  within  the  oonputational  grid  have  boundary 
conditicn  types  (BCIYPE)  associated  with  them. 

(5)  CALL  COT  =>  gives  COT  boundary  condition  input 
information  for  block-to-block  conditions  only.  MISSE  gets  back 
a  check  on  proper  boundary  alignment  and  phantom  point  locations. 

(6)  CALL  SEJIMEM  =>  gives  the  current  location  of  the  last 
memory  location  in  the  ribbon  vector  and  allocates  the  proper 
amount  of  memory.  [X-MP  version  only] 

(7)  CALL  POINT  =>  calculates  the  amount  of  memory  recjuired 
in  the  ribbon  vector. 

(8)  CALL  GRIDIN  =>  reads  the  computational  grid  from 
filetape  10  and  stores  it  in  the  parcper  memory  location  in  the 
ribbon  vector. 

(9)  CALL  IC  =>  calculates  the  initial  conditions  for  the 
d^jendent  variables  and  stores  them  in  the  proper  memory 
locations  in  the  ribbon  vector. 

(10)  CALL  READ  =  >  READ  will  read  the  grid  and  dependent 
variables  from  a  restart  file  only  if  one  is  available  and  is 
initiated  from  namelist  input  RESIKT  >  0  [X-MP  version  only],  or 
the  specified  restart  file  for  the  Cray  2  version. 

(11)  CALL  POTBC  =>  gives  PUIBC  the  memory  locations  of  the 
dependent  variables  and  then  puts  block-to-block  boundary 
condition  infornetion  in  manory  for  use  in  the  main  iteration 
loop. 

(12)  CALL  SEJIMPB  =>  fron  thte  boundary  condition  input,  SEIMPB 
sets  the  number  of  phantom  points  used  for  block-to-block 
calculations  and  stores  them  in  the  proper  memory  locations  in 
the  ribbon  vector. 

(13)  CALL  DPMAP  =>  sets  the  arrays  to  obtain  the  diagonal 
indices  of  the  matrices  for  vectorization. 

(14)  CALL  MEnuC  =>  calculates  the  areas  of  each  cell  face 
throughout  the  ccmputaticml  grid. 

(15)  CALL  OPSSD  =>  opens  an  SSD  file  and  assigns  it  the 
requested  number.  [X-MP  version  only] 

(16)  CALL  WRITE  =>  writes  file  information  to  SSD  [X-MP 
version  only]  or  to  the  specified  restart  file  [Cray  2  version]. 

Within  the  main  iteration  loop,  including  an  inner  loop  over  all 
the  blocks,  the  following  subroutines  are  called: 

( 1 )  CALL  POINT  =>  calculates  the  amount  of  memory  required  in 
the  ribbon  vector. 

(2)  CALL  SEIMEM  =>  sets  the  amount  of  memory  required  at  this 
point  in  the  iteration  locp  if  there  is  not  enough  allocated. 
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[X-MP  version  only] 

(3)  CAJjL  OPSSD  =>  opens  a  specific  SSD  file.  [X-MP  version 
only] 

(4)  CALL  READ  =>  reads  infornation  from  the  SSD  or  restart 
file  opaied  previously. 

(5a)  CALL  OPSSD  =>  opens  an  SSD  file  and  cissigns  it  a  specific 
nuttier  to  be  used  in  storing  the  flux  Jacobian  natrices  from  the 
positive  direction.  [X^IP  version  only] 

(5b)  CALL  OPSSD  =>  opens  an  SSD  file  and  cissigns  it  a  specific 
nuttber  to  be  used  in  storing  the  flux  Jacobian  natrices  from  the 
negative  direction.  [X-MP  version  only] 

(6)  CALL  SEIMEM  =>  sets  the  amount  of  memory  required  at  this 
point  in  the  iteration  loop  if  there  is  not  enough  allocated. 
[X-MP  version  only] 

(7)  CALL  READ  =>  reads  image  point  information  from  the  SSD  or 
restart  files. 

(8)  CALL  PVAR  =>  calculates  forces  and  mcments  for  each 
individual  surface  and  prepares  them  to  be  printed  out. 

(9)  CALL  BC  =>  sets  the  proper  boundary  conditions  for  each 
cell  on  the  boundaries. 

(10)  CALL  GBIBC  =>  obtains  the  boundary  condition  information 
for  block-to-block  Ccises  from  memory. 

(11)  CALL  SEIMEM  *>  sets  the  amount  of  memory  required  by  this 
portion  of  the  iteration  loc^  if  there  is  not  enough  allocated. 
[X-MP  version  only] 

(12)  CALL  STEP  =>  advances  the  flow  solution  one  iterative 
cycle. 

(13)  CALL  SEJMEM  =>  sets  the  amount  of  memory  required  by  this 
portion  of  the  iteration  loop  if  there  is  not  enou^  already 
allocated.  [X-MP  version  only] 

(14)  CALL  READ  =>  reads  image  point  information  from  the  SSD 
files. 

(15)  CALL  FUTBC  =>  puts  newly  calculated  block-to-block 
boundary  condition  information  in  memory  for  use  in  the 
appropriate  block. 

(16)  CauL  =>  writos  cut  all  information  in  the  ribbon 

vector  to  SSD.  [X-MP  version  only] 

At  the  end  of  the  main  iteration  loop,  the  Flow  Solver  then 
prints  out  the  forces  and  moments  calculated  during  the  entire 
run. 
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2.  SUBROtrriNE  ABQLU 


ABQLU  is  called  by  DOOM  (Section  4,7)  and  DOOP  (Section  4.8) 
to  corrpute  an  LU  deconposition  for  the  solution  matrix  on  the 
left  hand  side  of  Equation  (20).  AEQLU  receives  the  number 
of  necessary  points  in  the  grid  and  accotiplishes  an  LU 
deconposition  of  the  5X5  bloc)c  matrices.  ABQLU  returns 
the  deccnposed  matrices  through  COMMON  /E/. 

3.  SUBROUTINE  BC 

BC  is  called  by  the  main  routine,  MISSE,  to  apply  the 
appropriate  boundary  conditions  to  the  correct  patch  of 
corputational  grid.  BC  receives  information  from  the  Q 
array  (the  d^aendent  variables  =>  R,  RU,  RV,  RW,  E  and  pressure 
=  >  P),  the  metrics  (A,  SA),  and  the  size  of  the  blodc  (NI,  NJ, 

NK,  Nil,  NJ1,  NK1,  IB).  BC  returns  the  prcperly 
iirplemented  boundary  conditions  for  inpermeable  wall  (BCTyPE=1 ) 
and  farfield  conditions  (BCTyPE=2)  by  calling  SOLID  and 
FARFLD. 


4.  SUBROUTINE  CHECK 

CHECK  is  called  by  MISSE  from  within  the  read  BCIN  locp  and 
checks  each  line  of  boundary  condition  input  to  ensure  that  one 
and  only  one  boundary  condition  type  (BCTYPE)  is  specified  for 
each  corputational  cell.  CHECK  receives  all  of  the  BCIN 
input  information  (BLKA,  STARTA,  ENDA,  BLKB,  STARTB,  ENDB, 
BCIYPE),  the  maximum  values,  over  all  the  blocks,  for  NI  (NIX), 
NJ  (NJX),  and  NK  (NKX),  the  storage  locations  for  the  checking 
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planes  (IPL,  JPL,  KPL),  the  total  number  of  blocks  in  the  region 
(NBB).  CHECK  returns  the  internal  counter  for  the  number 
of  boundary  condition  cards  (NUMBC),  and  the  error  flag  (NERR  => 
set  equal  to  one  if  a  fatal  boundary  condition  error  is 
encountered). 

5.  SUBROUTINE  CUT 

CUT  is  also  called  by  the  main  routine  during  the  reading  of 
the  BCIN  inputs,  but  is  only  implemented  when  block-to-block 
conditions  are  raised.  The  CUT  subroutine  is  used  to  specify  the 
number  of  phantom  points  to  be  used  in  block-to-block 
calculations.  CUT  receives  the  BCIN  information  (BLKA, 

STARTA,  BNDA,  BLKB,  STARTB,  ENDB,  BCTYPE).  CUT  returns  the 
Image  point  information  (IMGE),  and  the  error  flag  (NERR  =>  set 
to  one  when  a  fatal  error  is  encountered).  Subroutine  COT 
calls  PASTE  (Section  4.24)  which  puts  the  boundary 
information  along  the  interfacing  boundary  (BLKB,  STARTB,  ENEB). 

6.  SUBROUTINE  DELQ 

DELQ  is  called  by  STEP  (Section  4.38)  to  calculate  the 
change  in  the  dependent  variables  using  the  Steger  axxd  Wanning 
flux-vector  splitting  technique  discussed  earlier  in  Section  2.2 
(SPLIT=1  ).  DELQ  receives  the  dependent  variables  and 
pressure  (R,  RU,  RV,  RW,  E,  P),  the  metrics  (AIX,  AJX,  AKX,  AIY, 
AJY,  AKY,  AIZ,  AJZ,  AKZ,  SADAI,  SADAJ,  SADAK),  the  number  of 
phantom  points  used  on  block  boundari'^s  (NPii,  NT) 2,  NP21,  NP22, 
NP31,  NP32),  and  the  indices  for  the  computational  block 
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currently  in  menory  (NI,  NJ,  NK,  Nil,  NJ1 ,  NK1,  IB).  DEHjQ 
returns  the  residual  vector  for  the  current  interation  (C®, 

DRU,  KIV,  [»w,  DE).  Subroutine  DELQ  calls  FLUX  (Section 
4.12)  to  obtain  the  flux  through  a  given  set  of  cell  faces 
(constant  I,  J,  or  K). 

7.  SUBROUTINE  DOOM 

DOOM  is  called  by  STEP  (Section  4.38)  to  solve  for  the  lower 
triangular  system  of  equations  using  Doolittle's  method.  DOOM 
receives  the  change  in  the  dependent  variables  (DR,  DRU, 

DRV,  DRW,  DE)  from  Equation  (21a)  for  the  previous  iteration,  the 
time  step  size  (DT),  the  coefficients  of  the  flux  Jacobian 
matrices  (B4,  B5,  B6),  the  size  of  the  conputational  bloc)c  (NI,  • 
NJ,  NK,  Nil,  NJ1,  NK1,  IB),  the  [MAPI,  MAP2]  information,  and  the 
SSD  file  to  read  the  flux  Jacobian  matrices  information  from 
(JFILE).  DOOM  returns  the  updated  changes  in  the  dependent 
variables  (CR,  DRU,  KIV,  DRW,  DE)  for  the  left  liand  side  of 
Equation  (21b).  Subroutine  D(XM  calls  ABQLU  to  obtain  the 

LU  decomposition  of  the  matrices. 


8.  SUBROUTINE  DOOP 

DOOP  is  called  by  STEP  (Section  4.38)  to  solve  for  the  upper 
triangulcur  system  of  equations  using  Doolittle's  method.  DOOP 
receives  the  residual  vector  (DR,  DRU,  DRV,  DRW,  DE)  for 
the  left  hand  side  of  Equation  (21a),  the  tine  step  size  (DT), 
the  coefficients  of  the  flux  Jacobian  matric‘-*s  (Bl,  B2,  B3),  the 
size  of  the  computational  bloc)c  (NI,  NJ,  NK,  Nil,  NJ1,  NKl ,  IB), 
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the  diagonal  plane  information  for  full  vectorization  (MAPI, 
MAP2),  and  the  SSD  file  to  read  the  flux  Jacobian  matrices 
information  from  (JFILE).  DOOP  returns  the  updated  changes 
in  the  dependent  variables  (DR,  DRU,  DRV,  ra?W,  DE)  for  the  left 
hand  side  of  Equation  (21b).  Subroutine  [XX)P  calls  AB2LU 
to  obtain  the  LU  decorposition  of  the  matrices. 

9.  SUBROOTINE  DPMAP 

DPMAP  is  called  by  the  main  routine  to  set  the  diagonal 
planes  in  the  solution  matrix  for  vectorization.  DFMAP 
coordinates  the  indices  interd^jendencies  so  as  to  allow  full 
vectorization  of  the  loops.  DHiAP  receives  the  size  of  the 
block  (NI,  NJ,  NK)  and  returns  the  diagonal  plane 
information  for  full  vectorization  (MAPI,  MAP2,  and  NDV). 

10.  SUBROUTINE  FARFLD 

EARFLD  is  called  by  BC  (Section  4.3)  to  conpute  the  boundary 
condition  information  required  at  the  outer  boundaries  of  the 
conputational  region.  EARFLD  receives  the  dependent 
variables  and  pressure  (R,  RU,  RV,  RW,  E,  and  P),  the  metrics  (A, 
SA),  the  size  of  the  block  (NI,  NJ,  NK,  Nil,  NJ1 ,  NK1 ,  IB),  and 
the  nunter  of  boundary  condition  input  cards  (M)  to  cycle  through 
to  check  for  farfield  patches  of  grid  (BCrYPE=2).  FARFLD 
returns  the  dependent  variables  at  the  outer  or  farfield 
boundctry  locations.  The  fcurfield  boundary  condition 
implementation  is  discussed  in  furth‘-?r  detail  in  Section  2.2. 
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1 1 .  SUBROUTINE  FJMAT 


FJMAT  is  called  by  STEP  (Section  4.38),  when  required  by  the 
specified  number  of  updates  to  the  Jacobians  (JFRBQ),  to 
calculate  the  coefficient  matrix  of  the  change  in  the  dependent 
variables  (the  flux  Jacobians).  FJMAT  receives  values  for 
the  indices  of  locps  (II,  J1,  Kl,  SIC34,  IB,  JB,  KB),  the  metrics 
(AX,  AY,  AZ,  SADA),  the  dependent  variables  and  pressure  (R,  RU, 
RV,  RW,  E,  and  P)  and  the  size  of  the  blodc  (NI,  NJ,  NK,  Nil, 

NJ1,  NK1,  IBLCXZK).  FJMAT  returns  updated  values  for  the 
flux  Jacobians  (B). 

12.  SUBROUTINE  FLUX 

FLUX  is  called  by  DELQ  (Section  4.6)  to  calculate  the  fluxes 
through  all  of  the  constant  I,  J,  or  K  cell  faces.  FLUX 
calculates  the  eigenvalues  associated  with  the  right  (EVIR,  EV4R, 
EV5R)  and  left  (EVIL,  EV4L,  EV5L)  hand  side  of  a  cell  face  to 
determine  the  "upwind”  side  for  differencing  the  dependent 
variables.  FLUX  receives  a  single  value  that  determines 
the  length  of  the  loop  and  is  tlte  number  of  points  in  either  the 
I  (NI),  J  (NJ),  or  K  (NK)  direction  and  FLUX  receives  the 
d^»ndent  variables  extrapolated  to  the  left  and  right  cell  face 
through  COMMON  /A/  and  /B/.  FLUX  returns,  through  CXDtMON 
/C/,  the  change  in  the  d^jendent  variables  for  the  right  hand 
side  of  Equation  (20)  (XR,  XRU,  XRV,  XRW,  XE). 
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13.  SUBROUTINE  FORCE 


FORCE  is  called  by  PVT^  (Section  4.28)  to  sum  the  forces  and 
moments  for  each  of  the  surfaces  specified  by  SURF,  SREF,  and 
Surname  and  to  set  up  the  print  format  for  the  output  of  the 
pressure  coefficients.  FORCE  receives  pressure  (P),  cell 
face  areas  (A),  axial  coordinate  locations  (X),  the  size  of  a 
given  block  (Nil,  NJ1,  NKl,  IB),  and  the  number  of  boundary 
condition  cards  (M)  to  cycle  thrcxigh  to  obtain  the  correct 
patches  of  grid  associated  with  iirpermeable  wall  boundary 
conditions  (BCTYPE=1 ).  FORCE  returns  pressure  coefficients 
(OOEP),  and  axial  cell  center  locations  (XC).  Subroutine  FORCE 
calls  GEPCP  to  obtain  the  pressure  coefficients  along  the 
solid  body  boundaries. 

14.  SIJBROUTINE  GETBC 

GEIBC  is  called  by  BC  (Section  4.3)  to  obtain  the  location 
of  the  block-to-block  boundary  condition  information.  GETBC 
receives  the  Q  array  information  (the  dependent  variables 
=>  R,  RU,  RV,  RW,  E  and  pressure  =>  P),  the  QBC  array  information 
(the  dependent  variables  on  the  block  boundaries  =>  QBC),  the  DQ 
array  (the  residuals  on  the  block  boundaries  =>  DR,  DRU,  DRV, 

DRW,  DE),  the  image  or  phantom  point  information  (IMGE),  and  the 
size  of  the  grid  (Nil,  RJ1,  NKl,  IB).  GETBC  returns  the 
necessary  location  of  the  blor'k-to-block  conditions  on  which  the 
Q  array  information  must  be  implemented.  Subroutine  GETBC 
calls  GETBBC  ( Section  4.15)  which  p5.ssGS  the  def'endent 
variables  (Q)  and  the  residuals  (DQ)  back  for  use  in  GETBC. 
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15.  SUBROUTINE  GETBBC 


GEIBBC  is  called  by  GBTBC  (Section  4,14)  to  obtain  the 
actual  block-to-block  boundary  information.  GEJTBBC 
receives  the  Q  and  DQ  arrays  (dependent  variables  =>  R,  RU, 

RV,  RW,  E,  pressure  =>  P,  residuals  =>  DR,  DRU,  DRV,  DRW,  DE), 
the  size  of  the  block  (Nil,  NJ1,  NK1,  IB),  and  the  image  point 
information  (IMGE,  N,  M).  GETBBC  returns  the  proper 
dependent  variables  and  residuals  for  use  in  GETBC  (Section 
4.14). 

16.  SUBROtrriNE  GETCP 

GETCP  is  called  by  FORCE  (Section  4.13)  to  calculate  the 
pressures  and  pressure  coefficients  along  all  of  the  solid  body 
boundaries  (BCTYPE=1 ).  GETCP  receives  the  size  of  the 
grid  (Nil,  NJ1 ,  NK1 ,  IB),  and  the  number  of  bounday  condition 
cards  (M)  to  sort  through  to  find  the  impermeable  wall  (BCrYPE=1 ) 
boundary  condition  patches  of  grid.  GETCP  returns  the 
pressures  (P)  and  pressure  coefficients  (COEP). 

17.  SUBROUTINE  GRIDIN 

GRIDIN  is  called  by  the  main  routine  to  read  in  the 
conputational  grid  tliat  discretizes  the  domain  about  the 
configuration  of  interest  and  to  calculate  the  X,  Y,  and  Z 
coordinates  of  the  grid  cells.  GRIDIN  receives  the  block 
number  in  which  to  calculate  the  coordinates  (IB),  and 
returns  the  coordinates  of  the  grid  (X,  Y,  Z).  GRIDIN 
reads  the  grid  file  off  of  filetape  10  (ETIO)  [X-MP  version 
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only]  and  file  1  for  the  Cray  2  version  of  the  code. 


18.  SUBROUTINE  IC 

IC  is  called  by  the  main  routine  to  obtain  the  initial 
conditiais  in  each  cell  throughout  the  entire  ccnputational 
domain.  IC  returns  the  dependent  variables  (R,  RU,  RV,  RW, 

E,  cind  P),  and  receives  the  size  of  the  grid  block  (Nil, 

NJl,  NK1,  IB).  The  dependent  variables  throu^K>ut  the  field  are 
initialized  bo  the  freestream  conditions  (RINIT,  RUINTT,  RVINTT, 
RWINIT,  EINIT,  and  PINTT)  which  are  received  through  OOMMCK 
/INIT/  and  CCMIOSI/CONST/ . 

19.  SUBROUTINE  LOCALT 

LOCALT  is  called  by  STEP  (Section  4.38)  to  obtain  the 
maximum  allowable  time  step  size  (DT)  for  each  cell  in  the  field. 
LOCALT  receives  the  dependent  variables  and  pressure  (R, 

RU,  RV,  RW,  E,  and  P),  the  metrics  (AIX,  AJX,  AKX,  AIY,  AJY,  AKY, 
AIZ,  AJZ,  AKZ,  SADAI,  SADAJ,  SADAK),  and  the  size  of  the  block 
(NI,  NJ,  NK,  Nil,  NJl,  NK1,  IB).  LOCALT  returns  the  value 
of  the  time  step  (DT)  for  each  cell  in  the  field  (Local  time- 
stuping)  which  is  calculated  using  the  speed  at  which  the 
prcpogation  of  information  is  moving  (C)  and  the  Courant  number 
(CFL)  input  by  the  user. 


20.  SUBROUTINE  MAIM 

MAIM  is  called  by  ROTX  (Section  4.32),  ROTY  (Section  4.33), 
and/or  ROTZ  (Section  4.34)  to  obtain  the  current  rotation  matrix 
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to  transform  the  initial  velocities  (RUINIT,  RVINIT,  PWINTT). 
MAIM  receives  a  transformation  matrix  (MAT)  from  ROIX, 

ROTY,  or  ROIZ  and  the  current  rotation  matrix  (T(3,3))  from 
OCmK  /TRANS/  and  returns,  through  OMOI  /TRANS/,  an 
updated  rotation  matrix  (T(3,3))  for  the  calculation  of  the 
initial  velocities. 

21 .  SUBROOTINE  METRIC 

METRIC  is  called  by  the  main  routine  to  calculate  the  areais 
of  the  cell  faces  (metrics)  throughout  the  entire  ccnputational 
region.  METRIC  receives  the  coordinates  of  the  grid  (X,  Y, 

Z),  and  the  size  of  the  grid  (NI,  NJ,  NK).  METRIC  returns 
the  areas  of  the  cell  faces  (AIX,  AJX,  AKX,  AIY,  AJY,  AKY,  AIZ, 
AJZ,  AKZ),  and  the  square  root  of  the  sum  of  the  squares  of  the 
metrics  (SADAI,  SADAJ,  SADAK). 

22.  SUBROOTINE  MINMOD 

MINMOD  is  called  by  RESID  (Section  4.30)  when  the  LIMIT=1 
option  is  selected.  MINMOD  is  a  flux  limiter  and  is  only 
enployeJ  when  the  Roe  flux-difference  splitting  technique 
(SPLIT=2)  is  selected.  MINMOD  receives  the  start  and  end 
indices  for  the  limiter  loop  (NSF,  NSL)  and  the  value  for  (PSI). 
MINMOD  returns  the  flux  limiter  (SP)  throug).  COMMCK 
/SIGMA/. 
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23.  SUBROOTINE  OPSSD  [X-MP  version  only] 


OPSSD  is  called  by  the  main  routine  to  open  SSD  files  for 
data  storage.  OPSSD  receives  the  unit  nunter  to  be  opened 
(lU)  and  the  name  of  the  SSD  file,  SSDissd  (ISSD).  Subroutine 
OPSSD  calls  routine  ASSIGN  to  cpen  and  assign  the  requested 
SSD  files. 

24.  SUBROUTINE  PASTE 

PASTE  is  called  by  CUT  (Section  4.5)  to  place  the  proper 
alignment  of  block-to-block  conditions  for  use  in  GBTBC  (Section 
4.14)  and  FUIBC  (Section  4.26).  PASTE  receives  the 
boundary  condition  information  necessary  to  properly  align  the 
block-to-block  interfaces  (ITORBA,  ICA,  ^KA,  STARTA,  ENDA,  ICB, 
BLKB,  STARTB,  NSC,  DC),  and  the  error  flag  (NERR).  PASTE 
returns  the  appropriate  interface  information  through  tlie 
image  arrays  (IMAGA,  M,  IMAGE,  N)  for  use  in  GETBC  and  FUTBC. 

25.  SUBROUTINE  POINT 

POINT  is  called  by  the  main  routine  to  calculate  the  length 
of  the  ribbon  vector.  POINT  receives  the  total  size  of  the 
blocks  (IS3D,  IS3D1,  IS3DP1,  IB)  and  the  current  location  of  the 
pointer  (NEXLOC).  POINT  returns  the  updated  location  of 
the  pointer  by  calculating  the  amount  of  memory  required  by  the 
information  needed  in  in-core  manory  (NEXLOC)  and  the  location  of 
the  starting  point  of  arrays  in  the  ribbon  vector  through  COMMON 
/FNTRS/. 
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26.  SUBROUTINE  PUIBC 


PUIBC  is  called  by  the  main  routine  to  place  the  proper 
d^endent  variables  in  the  QBC  arrays  for  use  in  block-to-block 
calculations.  PUTBC  receives  the  dependent  variables 
through  the  Q  array  (Q)  and  the  size  of  the  grid  (Nil ,  NJl ,  NK1 , 
IB)  and  returns  the  proper  values  in  the  QBC  arrays  (QBC). 
Subroutine  PUTBC  calls  PUIBBC  (Section  4.27)  to  obtain  the 
correct  dependent  variable  values  to  place  in  the  QBC  arrays. 

27.  SUBROUTINE  PUIBBC 

PUIBBC  is  called  by  PUIBC  (Section  4.26)  to  obtain  the 
correct  dependent  variables  to  place  in  the  QBC  arrays.  PUIBBC 
receives  the  dependent  variables  (Q),  the  size  of  the  blodc 
(Nil ,  NJl ,  NK1 ,  IB),  and  the  number  of  boundary  condition  cards 
to  cycle  through  to  coordinate  the  block-to-block  patches  of  grid 
(M).  PUIBBC  returns  the  QBC  array  with  the  correct  values 
for  the  block-to-block  interface  calculations  (QBC). 

28.  SUBROUTINE  FVAR 

PVAR  is  called  by  the  main  routine  to  calculate  the  forces 
and  ncments  for  all  of  the  surfaces  specified  by  SURF,  SREF,  and 
Surname.  PVAR  receives  the  pressure  (P),  the  areas  (A), 
the  axial  location  (X),  the  cell  center  locations  along  the  "x" 
axis  (XC),  the  number  of  boundary  condition  cards  for  obtaining 
the  s^lid  body  boundary  patches  of  grid  (M),  the  coefficient  of 
pressures  (COEP),  and  the  size  of  the  block  (NI,  NJ,  NK,  Nil, 

NJl,  NK1,  IB).  FVAR  returns  the  forces  and  moments  for 
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each  surface  through  COMMOI  /SEACE/.  Subroutine  PVAR  calls 
FORCE  (Section  4.13)  to  ctotain  the  intermediate  forces  and 


moments  for  each  surface. 

29.  SUBROUTINE  READ 

READ  is  called  by  the  main  routine  to  read  information  from 
the  restart  files.  READ  receives  the  specific  variables 
in-core  memory  location  (A),  the  amount  of  memory  required  by 
that  variable  (N),  and  the  SSD  unit  number  (lU).  READ 
returns  the  variable  itself  into  the  proper  memojry 
location. 

30.  SUBROUTINE  RESID 

RESID  is  called  by  STEP  (Section  4.38)  to  employ  the  Roe 
averaged  flux-difference  splitting  technique  discussed  in  Section 
2.3.  RESID  receives  the  dependent  variables  and  pressure 
(R,  RU,  RV,  RW,  E,  P),  the  metrics  (AIX,  AJX,  AKX,  AIY,  AJY,  AKY, 
AIZ,  AJZ,  AKZ,  SADAI,  SADAJ,  SADAK),  the  flux  limiter  selection 
(LIMIT),  the  phantom  point  information  for  block-to-block 
calculations  (NP11,  NP12,  NP21,  NP22,  NP31,  NP32),  and  the  size 
of  the  blocJc  (NI,  NJ,  NK,  Nil,  NJ1 ,  NK1 ,  and  IBLOCK).  RESID 
returns  the  ipdated  residuals  (DR,  DRU,  DRV,  DRW,  DE)  for 
the  right  hand  side  of  Equation  (21 ).  Subroutine  RESID 
calls  RLVECS  (Section  4.31  )  to  obtain  the  eigenvalues  and 
their  corresponding  right  and  left  eigenvectors  to  calculate  the 
flux,  and  calls  MINMOD  (Section  4.22)  or  SUPBEE  (Section 
4.39)  when  a  flux  limiter  is  required  by  the  LIMIT  variable  input 
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by  the  user.  For  the  X-MP  version  of  the  cxxte,  the  eigenvalues 
and  eigenvectors  are  stored  in  B1,  B2,  and  B3  shcuring  the  memory 
location  with  the  flux  Jacctoian  matrices.  For  the  Cray  2  version 
of  the  code,  the  eigenvalues  and  eigenvectors  are  siitply  stored 
in  locations  identified  by  EVL,  EVR,  and  EV. 


31 .  SUBROLfTINE  RLVEJCS 

RLVBCS  is  called  by  RESID  (Section  4.30)  to  calculate  the 
eigenvalues  and  their  corresponding  right  and  left  eigenvectors 
associated  with  the  right  hand  side  of  Equation  (20).  RLVBCS 
receives  loop  indices  (II,  J1,  K1,  IB,  JB,  KB),  the  metrics 
(AX,  AY,  AZ,  SADA),  the  d^)endent  variables  and  pressure  (R,  RU, 
RV,  RW,  E  and  P),  and  the  size  of  the  blodc  (NI,  NJ,  NK,  Nil, 
NJ1,  NK1,  IBLOCK).  RLVBCS  returns  the  calculated 
eigenvectors  and  eigenvalues  (B1,  B2,  B3  or  EVL,  EVR,  EV). 

32,  SUBROUTINE  ROTX 

ROIX  is  called  by  the  main  routine  to  calculate  a  revised 
rotation  matr  x.  ROIX  receives  the  magnitude  of  the 
rotation  about  the  "x”  axis  in  radians  (X)  and  iitplements  that 
angle.  ROTX  returns  a  revised  rotation  matrix  (T(3,3)) 
through  the  COMMIT  /TRANS/.  Subroutine  ROTX  calls  MATM 
(Section  4.20)  to  perform  the  matrix  multiplication  between  the 
original  transformation  matrix  and  the  requested  rotation 
(matrix)  to  yield  the  revised  rotation  matrix. 
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33.  SUEROtrriNE  ROIY 


ROIY  is  called  by  the  main  routine  to  calculate  a  revised 
rotation  matrix.  ROTY  receives  the  magnitude  of  the 
rotation  about  the  "y"  axis  in  radians  (X)  and  implements  that 
angle.  ROIY  returns  a  revised  rotation  matrix  (T(3,3)) 
through  the  OCMCN  /IRANS/.  Subroutine  ROIY  calls  MAIM 
(Section  4.20)  to  perform  the  matrix  multiplication  between  the 
original  transformation  matrix  and  the  requested  rotation 
(matrix)  to  yield  the  revised  rotation  matrix. 

34.  SUBROUTINE  ROTZ 

ROTZ  is  called  by  the  main  routine  to  calculate  a  revised 
rotation  matrix.  ROE2  receives  the  magnitude  of  the 
rotation  about  the  "z"  axis  in  radians  (X)  and  inplements  that 
angle.  ROTZ  returns  a  revised  rotation  matrix  (T(3,3)) 
through  the  ODNffW  /TRANS/-  Subroutine  ROTZ  calls  MAIM 
(Section  4.20)  to  perform  the  matrix  multiplication  between  the 
original  transformation  matrix  and  the  requested  rotation 
(matrix)  to  yield  the  revised  rotation  matrix. 

35.  SUBE^DUITNE  SETMEM  [X-MP  vei'sion  only] 

SETMEM  is  called  by  the  main  routine  to  set  the  memory 
requirements  for  the  flow  solver.  SETMEM  receives  the 
amount  of  memory  required  (MNEED)  and  returns  the 
sufficient  amount  of  memory  through  COMW  /MEM/.  Subroutine 
SEJIMEM  calls  the  system  routine  MEMORY  to  obtain  the 
required  amount  of  memory. 
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36.  SUBROUTINE  SEIMre 


SETMF©  is  called  by  the  main  routine  in  a  loop  over  all  of 
the  boundary  condition  cards  to  set  the  amount  of  idiantcm  points 
required  for  block-to-block  calculations.  SEHUPB  receives 
the  dimensions  of  the  particular  patch  of  block  (M2,  M3),  the 
indices  for  the  loop  over  that  patch  (IS,  JS,  IE,  JE),  and  the 
nurrber  of  phantom  points  required  (ILAP)  set  by  the  user  irpat 
BORDER.  SEIMPB  returns  the  MPB  array  (NP1 1 ,  NP12,  NP21 , 

NP22,  NP31,  NP32)  with  the  f^antom  point  information. 

37.  SUBROUTINE  SOLID 

SOLID  is  called  by  BC  (Section  4.3)  to  obtain  boundary 
corriiticxis  for  impermeable  wall  (BCrYPE=1 )  patches  of  grid. 

SOLID  receives  the  dependent  variables  and  pressure  (R,  RU, 

RV,  RW,  E,  and  P),  the  metrics  (A,  SA),  the  size  of  the  block 
(NI,  NJ,  NK,  Nil,  NJ1 ,  NK1 ,  IB),  and  the  number  of  boundary 
condition  cards  to  cycle  through  (M).  SOLID  returns  the 
values  of  the  dependent  variables  on  the  solid  bocfy  boundaries 
(R,  RU,  RV,  RW,  E). 

38.  SUBROUTINE  STEP 

STEP  is  called  by  the  main  routine  to  update  the  dependent 
variables.  STEP  receives  the  dependent  variables  and 
pressure  (R,  RU,  RV,  RW,  E  and  P),  the  metrics  (AIX,  AJX,  AKX, 
AIY,  AJY,  AKY,  AIZ,  AJZ,  AK2,  SADAI,  SADAJ,  SADAK),  the  flux 
limiter  selection  (LIMIT),  the  phantom  point  information  (NP11, 
NP12,  NP21,  NP22,  NP31 ,  NP32),  the  diagonal  plane  information  for 
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full  vector izat ion  (MAPI,  MAP2,  NDV),  the  splitting  technique  to 
be  utilized  (SPLIT),  the  size  of  the  block  (NI,  NJ,  NK,  Nil,  NJ1 , 
NK1,  IB),  and  the  print  cycles  (NCTC).  STEP  returns  the 
updated  dependent  variables  (R,  RU,  RV,  RW,  E  and  P).  Subroutine 
STEP  calls  either  RESID  (Section  4.30)  or  DELQ  (Section 

4.6)  to  obtain  the  residual  vector,  LCXZALT  (Section  4.19)  to 
obtain  a  value  for  the  time  step  size  (DT),  FJMAT  (Section  4.11 ) 
to  obtain  the  flux  Jacobian  matrices  (A"*",  A”,  B'*', 

B~,  C*",  C“),  WRITE  (Section  4.40)  to  write  out  the 

flux  Jaoobians  to  SSD,  and  DOOP  (Section  4.8)  and  DOOM  (Section 

4.7)  to  solve  the  block  5X5  upper  or  lower  system  of  equations. 

39.  SUBROUTINE  SUPBEE 

SUPBEE  is  called  by  RESID  (Section  4.30)  to  obtain  flux 
limiters  when  the  LIMIT=2  option  is  selected  by  the  user.  SUPBEE 
receives  the  loop  indices  (NSF,  I^L)  and  returns  the 
flux  limiters  through  CDM01  /SIQiA/. 

40.  SUBROUTINE  WRITE 

WRITE  is  called  by  the  main  routine  to  write  information  out 
to  restart  files.  WRITE  receives  the  in-core  memory 
location  of  the  specified  variable  (A),  and  the  amount  of  memory 
required  by  the  Vctriable  (N).  WRITE  returns  the  specified 
variable  into  the  restart  SSD  file  (lU). 
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SBCnCN  V 


APPLICATICNS 

The  Program  EAGLE  -  Flow  Solver  has  been  extensively 
exercised  on  several  aerodynamic  problems  to  validate  the  codes 
predictive  capabilities.  These  cases  have  been  run  through  a 
range  of  Mach  numbers  (compressible  subsonic,  transonic, 
supersonic)  and  incidence  angles  (0  to  20  degrees).  The  Flow 
Solver  yields  excellent  inviscld  results  when  ccnpared  to 
experimental  data  and  should  provide  a  very  reasonable  estimate 
of  the  pressure  distribution,  as  well  as  forces  and  moments  for 
any  arbitrary  configuration.  Presented  in  this  section  are  two 
sanples  from  this  validation  process  that  will,  hopefully,  aid  in 
the  understanding  of  the  uses  of  the  Program  EAGLE  -  Flow  Solver. 
A  discussion  is  provided  for  each  grid,  describing  the  single  or 
multi-block  structure  and  the  cxarputational  region  on  which  the 
Flow  Solver  actually  works.  The  Flow  Solver  input  for  each 
example  case  is  presented  with  a  discussion  of  the  individual 
input  statements. 

1 .  NACA  0012  AIRFOIL 

The  first  sanple  application  involves  a  syirmetric  airfoil 
(NACA  0012)  at  transonic  conditions  (Mach  =  0.80)  at  a  small 
incidence  (1.25  degrees).  This  case  has  been  a  primary  check 
case  throughout  the  development  of  the  EAGLE  code  and  has  yielded 
some  interesting  results  when  conparing  the  different  splitting 
techniques  ( Steger-Warming  =>  SPLIT=1  and  Roe  =>  SPLIT=2). 
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a.  One  Block  Grid 

Ihe  Program  EAGLE  -  Numerical  Grid  Generation  System  has 
been  enployed  to  generate  a  single  block  algebraic,  two- 
dimensional  C-type  grid  about  the  NACA  0012  airfoil  (Figure  10a). 
The  221  X  20  X  2,  C-type  grid  is  used  so  as  to  allow  a  fine 
concentration  of  grid  points  around  the  leading  edge  of  the 
section  (Figure  10b).  Grid  lines  have  been  concentrated  parallel 
(J  lines)  to  the  airfoil  body,  and  forward  near  the  leading  edge 
perpendicular  (I  lines)  to  the  airfoil  body,  to  resolve  the  sharp 
gradients  in  the  flowfield  expected  in  these  regions.  A  view  of 
the  corputational  region  is  shown  (Figure  11 ),  with  dimensions, 
to  help  with  the  input  of  boundary  conditions  (BCIN). 

b.  Flow  Solver  Input 

The  Program  EAGLE  -  Flow  Solver  input  is  designed  to  be  user 
oriented  in  that  each  namelist  iiput  has  a  specific,  readable, 
meaning  for  the  user.  In  this  case  the  input  stream  is  a  rather 
slitple  one  in  that  the  two-dimensional  grid  requires  little 
boundary  oordLLtion  input.  The  input  for  this  case  follows: 


E$  FINPUT 

FSMACH  =  0.80,  CFL  =  15.0  , 

GAM  =  1.4 

f 

NGRAD  =  0  ,  NZPG  =  10  , 

JFRBQ  =  20 

f 

NPR  =2  ,  NIT  =  500  , 

PRINT  =  1 

f 

BORDER  =  2  ,  SPLIT=  2 

LIMIT  =  2 

t 

LIFTAX  =2  $ 

E$  ROTATE 

TRANS  =  3  ,  PHIZ  =  -1 .25 

$ 

NACA  0012 

E$  SURFACE  SREF  =  1  ,  XREF  =  0.25 

$ 

The  three  general  input  types  are  presented  here  for  ease  of 
discussion.  In  the  FINPUT  section;  the  freestream  Mach  nunfcer  is 
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set  to  0.80  (FSMACH),  the  Courant-Friedrich-Lewy  (CFL)  cxmdition 
is  set  at  15.0,  the  ratio  of  specific  heats  (GAM)  1.4  for  air. 
There  will  be  10  zero  pressure  gradient  boundary  conditions 
(NZPG)  and  none  will  be  gradually  applied  (NGRAD  =0).  The  flux 
Jacobian  matrices  will  be  updated  every  20  iterations  (JFT^)  and 
two  phantom  points  will  be  used  between  blocks  for  block-to-block 
conditions  (BORDER).  The  total  number  of  iterations  will  be  1000 
(NPR*NIT)  and  there  will  be  five  printouts  (NPR)  of  pressures 
versus  X  location,  forces  and  moments  (PRINT)  every  500 
iterations  (NIT).  The  splitting  technique  ertployed  (SPLIT)  will 

be  a  form  of  Roe^  averaged  flux-difference  splitting  (Section 
2.3)  using  the  SUPBEE  (LIMIT)  routine  (Section  4.39).  And 
finally,  the  direction  of  the  lift  axis  is  set  to  tiie  positive 
"y"  axis  (LIFTAX)  for  use  with  the  ROTATE  input  card. 

The  ROTATE  input  is  straightforward  in  that  it  describes  a 
rotation  about  the  positive  "z"  axis.  The  ”z"  axis  is  ret  as  the 
axis  of  rotation  (TRANS)  and  the  magnitude  of  this  angle  is 
described  as  -1.25  degree  (PHIZ).  Figure  12  shows  the  details 
of  this  rotation.  Since  all  rotations  are  counter-clodcwise 
positive,  and  tail  down  is  a  positive  angle  of  attack,  to  obtain 
a  1.25  degree  angle  of  attack  the  rotations  must  be  as  indicated. 

The  first  SURFACE  input  is  a  string  of  a  maximum  of 
sixteen  alphanumeric  characters  (NACA  0012).  This  aids  in  the 
description  of  the  surfaces  to  be  specified.  Following  this  card 
is  the  SURFACE  input  itself  in  which  most  of  the  default  values 
are  exercised  ( LREF  =  1.0  for  chord  length ) .  The  surface 
specification  throughout  this  run  will  be  1  (SREF)  for  the  entire 


airfoil  and  will  coincide  with  values  of  1  for  the  SURF  input  on 
the  BCIN  cards  (BCTYPE=1 )•  Ihe  reference  "x"  location  is  set  to 
the  quarter  chord  (XREF)  for  force  and  ncinent  calculations. 

The  following  cards  detail  the  BCIN  input  to  the  Flow  Solver: 

E$  BCIN  BCTYPE  =  1,  SURF  =  1 ,BLKA  =  1,  STARTA  =  31,1,  1, 

ENDA=  191,  1,  2  $ 

E$  BCIN  BCTYPE  =  2,  BLKA  =  1 ,  STARTA  =  1 ,  1 ,  1 , 

ENDA  =1,  20,  2  $ 

E$  BCIN  BCTYPE  =  2,  BLKA  =  1,  STARTA  =  221,1,1, 

ENDA  =  221,20,  2  $ 

E$  BCIN  BCTYPE  =  2,  BLKA  =  1,  STARTA  =  1,  20,1, 

ENDA  =  221,20,  2  $ 

E$  BCIN  BCTYPE  =  4,  BLKA  =  1 ,  STARTA  =  1 ,  1 ,  1 , 

ENDA  =31,1,2, 

BLKB  =  1,  STARTS  =  221,1,1, 

ENT®  =  191,1,  2  $ 

E$  BCIN  BCTYPE  =  0  $ 

For  this  exanple,  inpermeable  surface  (solid  body)  boundary 
conditions  (BCTYPE)  are  applied  in  blocJc  1  (BLKA)  from  31,1,1  to 

191.1.1  (STARTA  and  ENDA)  a  constant  J  line.  This  portion  of 
solid  body  will  be  used  in  the  force  and  moment  calculations  for 
surface  1  (SURF  =>  SREF  and  NACA  0012).  Far field  boundary 
conditions  (BCTYPE)  are  applied  in  block  1  (BLKA)  for  three 
different  planes  (lines,  in  this  two-dimensional  exanple)  from 

1.1.1  to  1,20,2,  a  constant  I  line  (the  upper  rear  boundary); 
from  221,1,1  to  221,20,2,  another  constant  I  line  (the  lower  rear 
boundary)  and  from  1,20,1  to  221,20,2,  a  constant  J  line  (the 
outer  boundary).  Block-to-block  conditions  are  applied  in  the 
wake  region  of  this  airfoil  in  block  1  (BLKA  and  BLKB)  from  1,1,1 
to  31,1,2  (STARTA  and  ENDA)  corresponding  to  a  patch  (line  in  2- 
D )  of  9-  -cd  from  221,1,1  to  191,1,2  ( STARTB  and  ENDB ) ;  again ,  a 
constant  J  line.  Finally,  the  BCIN  input  dec':  is  terminated  by 
employing  the  ending  inpxit  of  zero  (BCTYPE). 
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c.  Flow  Solver  Output 

Oul^t  for  the  Program  EftGLE  -  Flow  Solver  is  organized  to 
echo  all  user  inputs  for  checking  purposes,  print  convergence 
information  (residuals  and  number  of  supersonic  points),  print 
pressure  coefficient  data  (if  requested)  and  to  print  forces, 
moments  and  reference  lengths  and  areas.  Table  1  shows  the  output 
for  this  test  case.  Echo  checking,  from  the  namelist  input 
streams  and  from  specific  information  regarding  solution 
procedures,  is  displayed  as  shown.  The  dimensions  of  the 
numerical  grid  generated  are  echoed  for  checking  purposes  and  a 
header  is  displayed  giving  the  type  of  splitting  performed-  The 
convergoTce  information  by  iteration  number  (STEP)  and  block 
number  (BLCXZK)  is  shown.  Residuals  (RMS  and  MAX)  are  given  for 
density  and  energy. 

The  Flow  Solver  output  for  this  sanple  case  is  given  in  Table 
1  and  a  plot  of  the  pressure  coefficient  ((^)  versus  chord  length 
(X/C)  is  shown  in  Figure  16.  The  convergence  for  this  example  is 
quite  good  in  that  the  L2  Norm  (not  shown)  has  converged  to  -2.88 
or  three  orders  of  magnitude  in  1000  iterations,  the  residuals 
for  density  are  at  E-4  and  the  number  of  supersonic  points  is 
constant  at  451  for  twelve  iterations  (although  not  shown,  NSUP 
was  constant  at  450  for  almost  200  iterations).  All  indications 
show  that  a  good  engineering  level  solution  has  been  attained  in 
1000  iterations.  Any  further  convergence  would  go  essentially 
unnoticed  in  any  graphical  analysis  of  the  solution.  The  maximum 
memory  size  was  under  1.7  million  words  (maximum  field  length) 
and  the  run  time  was  115.335  CP  seconds  on  the  Cray  X-MP. 


TABLE  I.  Program  EAGLE  -  Flow  Solver 
NAO^  0012  Airfoil  Output 


*  * 

*  Program  EAGLE-  Flow  Solver  * 

*  Multi-block  Inplicit  Steady-State  Euler  Algorithm* 

*  * 

*  Cray  2  version  (20  July  1988)  * 

*  * 
************<r**************************ik************ 

Air  Force  Armament  Laboratory  (AFATL) 
Aeromechanics  Division  ( AFATL /FX) 

Aerodynamics  Branch  (AEATL/EXA) 

Gcnputational  Fluid  Dynamics  Section  (CFD) 

Eglin  AFB,  FL  32542-5242 

Questions  or  Coitinents  call: 

Jon  S.  Mounts  or  Dave  M.  Belk  @ 

(904)  882-3124/2767  or  AV  872-3124/2767 


******************  echo  check  ********************** 

Namelist  /FINPOT/  : 

Courant  number  (CFL)  =  15.0 
Freestream  Mach  number  (FSMACH)  =  0.8 
Frequency  of  flux  Jacobian  matrix  updates 
(JFRBQ)  =  20 

Number  of  zero  pressure  gradient  BC  iterations 
(NZPG)  =  10 

Number  of  gradually  applied  ZPG  BC  iterations 
(NGRAD)  =  0 

Lift  will  be  calculated  for  the  #2  axis 
(LIFIAX) 

Total  number  of  surfaces  (NSURF)  =  1 

Total  number  of  iterations  (NFR)*(N1T)  =  1000 

*cW*^t€W***************************************** 

Computational  region  : 

Block  NI  NJ  NK 

1  221  20  1 

Block  1  assumed  to  be  two-dimensional 

*************************************************** 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 

IticititiciilkitititHititltleititltlc’kiiic'kititicitltititlt'kiifkitifk'k'kit'k-kitic'kicieicitic 

Namelist  /ROTATE/  : 

Rotation  angle  sequence  (TRANS)  =3,0,0 
Rotation  angles  : 

PHIX  =  0. 

PHIY  =  0. 

PHIZ  =  -1.25 

************************************************* 

Differencing  Technique  : 

SPLIT  =  2 
LIMIT  =  2 
RORDER=  2 

************************************************** 

Namelist  /SURFACE/  : 

Surface  reference  number  (SREF)  :  1 

Reference  area  (AREF)  =  1.0 

Reference  length  (LREF)  =  1.0 

"x”  reference  (XREF)  =  0.25 

"y"  reference  (YREF)  =  0. 

"z"  reference  (ZREF)  =  0. 

ififkitltititititititititltlticititicititltltit'kititititit'kiiiiititicitititieicitititlt'kititic 

Total  memory  allocated  :  1664816 


2  PRINTOUIS  EVERY  500  CYCLES 


STEP 

IB  RTMAX 

IRM  JE1M  KRM  RIRMS 

EIMAX 

EIRMS 

0 

1  6.420e-01 

111 

3 

2  0. 

1 .284e^-00 

0. 

1 

1  3.819e-01 

111 

2 

2  -2.323e-01 

8.227e-01  -2.309e-01 

2 

1  1.752e-01 

113 

4 

2  -4.376e-01 

3.581e-01  -4.619e-01 

3 

1  -1.784e-01 

95 

2 

2  -4.249e-01 

-3.565e-01  -4.773e-01 

4 

1  -2.303e-01 

98 

2 

2  -3.450e-01 

-5.230e-01  -3.825e-01 

5 

1  -2.706e-01 

111 

2 

2  -3.133e-01 

-6.766e-01  -3.346e-01 

6 

1  -2.649e-01 

111 

2 

2  -3.354e-01 

-6.775e-01  -3.425e-01 

7 

1  -2.094e-01 

112 

4 

2  -4.051e-01 

-4.487e-01  -3.983e-01 

8 

1  -1.617e-01 

112 

5 

2  -4.733e-01 

-3.672e-01  -4.688e-01 

9 

1  1.405e-01 

97 

2 

2  -4.960e-01 

-3.263e-01  -5.113e-01 

10 

1  3.739e-01 

111 

2 

2  -2.892e-01 

8.077e-01  -3.253e-01 

11 

1  -1.327e-01 

95 

2 

2  -5.449e-01 

-2.983e-01  -5.505e-01 

12 

1  -1.158e-01 

96 

2 

2  -6.118e-01 

-2.549e-01  -6.408e-01 

13 

1  9.642e-02 

89 

5 

2  -6.515e-01 

-2.016e-01  -6.999e-01 

14 

1  1.006e-01 

88 

6 

2  -6.781e-01 

1.772e-C1  -7.330e-01 

15 

1  9.985e-02 

90 

6 

2  -6.993e-01 

1.854e-01  -7.516e-01 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 


485 

1  -1.028e-03 

60 

2 

2 

486 

1  -1.008e-03 

60 

2 

2 

487 

1  1.016e-03 

197 

2 

2 

488 

1  -1.001e-03 

60 

6 

2 

489 

1  -1.022e-03 

60 

6 

2 

490 

1  -9.777e-04 

60 

2 

2 

491 

1  -9.636e-04 

60 

2 

2 

492 

1  1.167e-03 

196 

2 

2 

493 

1  -9.166e-04 

60 

2 

2 

494 

1  -9.001e-04 

60 

6 

2 

495 

1  -8.905e-04 

60 

6 

2 

496 

1  1.060e-03 

195 

2 

2 

497 

1  9.354e-04 

195 

2 

2 

498 

1  -8.659e-04 

60 

2 

2 

499 

1  -8.445e-04 

60 

2 

2 

ISTEP 

=  500 

BLOCK  =  1 
2 

X 

CP 

32 

0.9948 

0.3623 

33 

0.9842 

0.3165 

34 

0.9729 

0.2636 

35 

0.961 1 

0.2231 

36 

0.9487 

0.1907 

37 

0.9357 

0.1620 

38 

0.9220 

0.1366 

39 

0.9078 

0.1130 

40 

0.8929 

0.0917 

41 

0.8773 

0.0711 

42 

0.8612 

0.0515 

43 

0.8445 

0.0320 

44 

0.8271 

0.0133 

45 

0.8092 

-0.0054 

46 

0.7907 

-0.0241 

47 

0.7717 

-0.0421 

48 

0.7521 

-0.0591 

49 

0.7321 

-0.0791 

50 

0.7117 

-0.1008 

51 

0.6909 

-0.1145 

52 

0.6698 

-0.1244 

53 

0.6483 

-0.1503 

54 

0.6266 

-0.1821 

55 

0.6048 

-0.1852 

56 

0.5828 

-0.1820 

57 

0.5608 

-0.2137 

58 

0.5387 

-0.2838 

59 

0.5168 

-0.3219 

60 

0.4949 

-0.3050 

61 

0.4732 

-0.2754 

62 

0.4517 

-0.2841 

63 

0.4305 

-0.3523 

64 

0.4097 

-0.3981 

65 

0.3892 

-0.4125 

66 

0.3692 

-0.4066 

67 

0.3496 

-0.3351 

68 

0.3305 

-0.6208 

69 

0.3120 

-0.6733 

-2.933e+00  -1.939e-03  -3.015ef00 
-2.943e+00  -1.928e-03  -3.025e+00 
-2.939e+00  -l-905e-03  -3.021e4-00 
-2.946e+00  -1.935e-03  -3.026e+00 
-2.949e+00  -1.971e-03  -3.029ef00 
-2.953ef00  -1.909e-03  -3.034e+00 
-2.949ef00  -1.867e-03  -3.033e+00 
-2.942ef00  1.916e-03  -3.029e+00 
-2.959ef00  -1.753e-03  -3.045e+00 
-2.968e+00  -1.700e-03  -3.053e+00 
-2.970e+00  -1.666e-03  -3.054ei-00 
-2.964e4-00  1.791e-03  -3.049ef00 
-2.972e+00  -1.628e-03  -3.054e+00 
-2.980e+00  -1.628e-03  -3.061e+00 
-2-977e4-00  -1  .574e-03  -3.060e+00 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 


70 

0.2940 

-0.6638 

71 

0.2766 

-0.6599 

72 

0.2598 

-0.6523 

73 

0.2437 

-0.6421 

74 

0.2281 

-0.6291 

75 

0.2132 

-0.6141 

76 

0.1989 

-0.5980 

77 

0.1852 

-0.5804 

78 

0.1722 

-0.5617 

79 

0.1598 

-0.5416 

80 

0.1480 

-0.5197 

81 

0.1368 

-0.4958 

82 

0.1262 

-0.4703 

83 

0.1161 

-0.4456 

84 

0.1066 

-0.4191 

85 

0.0977 

-0.3919 

86 

0.0893 

-0.3608 

87 

0.0813 

-0.3324 

88 

0.0739 

-0.3004 

89 

0.0669 

-0.2650 

90 

0.0603 

-0.2318 

91 

0.0542 

-0.1942 

92 

0.0485 

-0.1523 

93 

0.0432 

-0.1144 

94 

0.0382 

-0.0648 

95 

0.0336 

-0.0192 

96 

0.0294 

0.0348 

97 

0.0254 

0.0911 

98 

0.0218 

0.1534 

99 

0.0185 

0.2202 

100 

0.0155 

0.2929 

101 

0.0127 

0.3712 

102 

0.0103 

0.4598 

103 

0.0081 

0.5546 

104 

0.0062 

0.6547 

105 

0.0045 

0.7572 

106 

0.0031 

0.8591 

107 

0.0020 

0.9495 

108 

0.0012 

1.0346 

109 

0.0006 

1.0891 

110 

0.0002 

1.1097 

111 

0.0000 

1.1033 

112 

0.0000 

1.0684 

113 

0.0002 

1.0021 

1 14 

0.0006 

0.9023 

115 

0.0012 

0.7866 

116 

0.0020 

0.6513 

1 17 

0.0031 

0.5147 

1 18 

0.0045 

0.3785 

119 

0.0062 

0.2564 

120 

0.0081 

0.1422 

121 

0.0103 

0.0380 

122 

0.0127 

-0.0546 

123 

0.0155 

-0.1366 

124 

0.0185 

-0.2079 

125 

0.0218 

-0.2726 

126 

0.0254 

-0.3395 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 

127  0.0294  -0.3805 

128  0.0336  -0.4327 

129  0.0382  -0.4860 

130  0.0432  -0.5126 

131  0.0485  -0.5503 

132  0.0542  -0.5936 

133  0.0603  -0.6298 

134  0.0669  -0.6575 

135  0.0739  -0.6821 

136  0.0813  -0.7087 

137  0.0893  -0.7379 

138  0.0977  -0.7655 

139  0.1066  -0.7890 

140  0.1161  -0.8090 

141  0.1262  -0.8291 

142  0.1368  -0.8522 

143  0.1480  -0.8756 

144  0.1598  -0.8971 

145  0.1722  -0.9157 

146  0.1852  -0.9318 

147  0.1989  -0.9471 

148  0.2132  -0.9641 

149  0.2281  -0.9818 

150  0.2437  -0.9989 

151  0.2598  -1.0146 

152  0.2766  -1.0287 

153  0.2940  -1.0413 

154  0.3120  -1.0527 

155  0.3305  -1.0635 

156  0.3496  -1.0743 

157  0.3692  -1.0855 

158  0.3892  -1.0963 

159  0.4097  -1.1063 

160  0.4305  -1.1152 

161  0.4517  -1.1228 

162  0.4732  -1.1291 

163  0.4949  -1.1345 

164  0.5168  -1.1391 

165  0.5387  -1.1433 

166  0.5608  -1.1472 

167  0.5828  -1.1489 

168  0.6048  -1.1284 

169  0.6266  -0.1077 

170  0.6483  0.0132 

171  0.6698  0.0127 

172  0.6909  0.0122 

173  0.7117  0.0163 

174  0.7321  0.0208 

175  0.7521  0.0265 

176  0.7717  0.0345 

177  0.7907  0.0447 

178  0.8092  0.0563 

179  0.8271  0.0694 

180  0.8445  0.0833 

181  0.8612  0.0984 

182  0.8773  0.1145 

183  0.8929  0.1316 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 

184  0.9078  0.1509 

185  0.9220  0.1721 

186  0.9357  0.1960 

187  0.9487  0.2232 

188  0.9611  0.2552 

189  0.9729  0.2939 

190  0.9842  0.3422 

191  0.9948  0.3704 

500  1  9.415e-04  194  2  2  -2.979ef00  1.637e-03  -3.062ef00 

MACH  NO.  PHIX  PHIY  PHIZ  CFL  LIFIAX 

0.8000  0.00  0.00  -1.25  15.0  2 

REFERENCE 

SURFACE  NAME  AREA  X  Y  Z 

1  NACA  0012  Test  1.00000  0.25000  0.00000  0.00000 


SURFACE  NAME  FORCES 

1  NACA  0012  Test  ETC  = 

FY  = 

FZ  = 


SURFACE  NAME  FORCES 

1  NACA  0012  Test  CX  = 

a  = 

CZ  = 

SURFACE  NAME  FORCES 

1  NACA  0012  Ttest  FL  = 

FD  = 

FS  = 


501 

1  1.018e-03 

194 

2 

2 

502 

1  7.637e-04 

57 

8 

2 

503 

1  8.444e-04 

57 

8 

2 

504 

1  1 . 045e-03 

197 

2 

2 

505 

1  -1.009e-03 

195 

2 

2 

506 

1  8.323e-04 

58 

5 

2 

507 

1  8.435e-04 

58 

5 

2 

508 

1  8.424e-04 

58 

5 

2 

509 

1  l.173e-03 

196 

2 

2 

510 

1  8.417e-04 

58 

6 

2 

51 1 

1  8.434e-04 

58 

6 

2 

512 

1  8.403e-04 

58 

6 

2 

513 

1  8.933e-04 

195 

2 

2 

514 

1  1.007e-03 

195 

2 

2 

515 

1  3.197e-04 

58 

7 

2 

516 

1  8.579e-04 

58 

7 

2 

517 

1  8.385e-04 

58 

7 

2 

518 

1  1 . 076e-03 

194 

2 

2 

519 

1  1 . 072e-03 

194 

2 

2 

520 

1  7.381e-04 

58 

6 

2 

521 

1  8.580e-04 

197 

2 

2 

522 

1  -1.048e-03 

195 

2 

2 

MOMENTS 

0.009043  MX  =  -0.163014 
0.326028  MY  =  0.004522 
0.000000  MZ  =  0.029244 

COEFFICIENTS 

MOMENTS 

0.009043  CMX  =  -0.163014 
0.326028  CMY  =  0.004522 
0.000000  CMZ  =  0.029244 

COEFFICIENTS 

0.325753  CL  =  0.325753 
0.016153  CD  =  0.016153 
0.000000  CS  =  0.000000 

-2.977e+00  1.744e-03  -3.060e+00 
-2.991e^-00  1.540e-03  -3.074e+00 
-2.995e+00  1.617e-03  -3.078e+00 
-2.992e+00  1.778e-03  -3.076ef00 
-3.003e+00  -1.813e-03  -3.087e+00 
-3.008e4-00  1.584e-03  -3.092e+-00 
-3.006e+00  1.608e-03  -3.094e4-00 
-3.004e+00  1.622e-03  -3.093ei-00 
-3.000e+00  1.968e-03  -3.090e+00 
-3.014e+00  1.585e-03  -3.103ef00 
-3.018e+-00  1  .572e-03  -3.107ef00 
-3.017e+00  1.567e-03  -3.107ef00 
-3.019e+00  1.566e-03  -3.110e+00 
-3.023e+00  1.725e-03  -3.115e+00 
-3.031e+00  1.610e-03  -3.120ei-00 
-3.027e+00  1.660e-03  -3.112e+00 
-3.022e+00  1 .648e-03  -3. 106e+00 
-3.011 e+00  1  - 942e-03  -3 . 093e+00 
-3.023e+00  1  - 905e-03  -3. 1 05e+00 
-3 . 037e+00  1.481 e-03  -3.11 8e+00 
-3.028e+00  1.452e-03  -3.111ef00 
-3. 031 e+00  -1.879e-03  -3.114e+00 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 


-03  195  2  2 


524 

1  7.364e-04 

58 

7 

2 

525 

1  7.258e-04 

58 

7 

2 

975 

1  -7.631e-04 

197 

2 

2 

976 

1  8.539e-04 

195 

2 

2 

977 

1  9.282e-04 

195 

2 

2 

978 

1  3.958e-04 

27 

3 

2 

979 

1  -7.310e-04 

196 

2 

2 

980 

1  7.057e-04 

194 

2 

2 

981 

1  9.820e-04 

194 

2 

2 

982 

1  8.784e-04 

194 

2 

2 

983 

1  7.618e-04 

197 

2 

2 

984 

1  -9.734e-04 

195 

2 

2 

985 

1  -8.540e-04 

195 

2 

2 

986 

1  -4.825e-04 

198 

2 

2 

987 

1  -8.122e-04 

194 

2 

2 

988 

1  -1.218e-03 

194 

2 

2 

989 

1  9.981e-04 

196 

2 

2 

990 

1  -6.028e-04 

193 

2 

2 

991 

1  -6.008e-04 

193 

2 

2 

992 

1  -7.360e-04 

197 

2 

2 

993 

1  9.029e-04 

195 

2 

2 

994 

1  8.924e-04 

195 

2 

2 

995 

1  -4.598e-04 

193 

2 

2 

996 

1  -7.557e-04 

196 

2 

2 

997 

1  -7.188e-04 

196 

2 

2 

998 

1  7.838e-04 

194 

2 

2 

999 

1  9.871e-04 

194 

2 

2 

ISTEP  =  1000  BLOCK  =  1 
2 


X 

CP 

32 

0.9948 

0.3675 

33 

0.9842 

0.3198 

34 

0.9729 

0.2663 

35 

0.961 1 

0.2258 

36 

0.9487 

0.1930 

37 

0.9357 

0.1643 

38 

0.9220 

0.1392 

39 

0.9078 

0.1153 

40 

0.8929 

0.0934 

41 

0.8773 

0.0724 

42 

0.8612 

0.0528 

43 

0.8445 

0.0337 

44 

0.8271 

0.0153 

45 

0.8092 

-0.0028 

46 

0.7907 

-0.0206 

47 

0.7717 

-0.0386 

48 

0.7521 

-0.0565 

49 

0.7321 

-0.0747 

50 

0.7117 

-0.0930 

51 

0.6909 

-0.1115 

-3.034e4-00  -1.917e-03  -3.115e+00 
-3.041ei-00  1.43le-03  -3.120e+00 
-3.044e+00  1-480e-03  -3.126e+00 


-3.314ef00  -1.345e-03  -3.420e+-00 
-3.339ei-00  1.434e-03  -3.446e+00 
-3.343e4-00  1.592e-03  -3.442ef00 
-3.352e+00  7.469e-04  -3.445ef00 
-3.320e+00  -1.168e-03  -3.419e+00 
-3.343ef00  1.119e-03  -3.436ef00 
-3.343e4-00  1.648e-03  -3.435e+00 
-3.345e+00  1.510e-03  -3.436e4-00 
-3.331e+00  1.27le-03  -3.436e4-00 
-3.319ei-00  -1.676e-03  -3.419e+00 
-3.348ef00  -1 .476e-03  -3.440e+00 
-3.390ef00  -9.262e-04  -3.496e+00 
-3.374ef00  -1.316e-03  -3.491e+00 
-3.297e4-00  -2.012e-03  -3.413e+00 
-3.313e4-00  1.70le-03  -3.418e+00 
-3.369e+00  -9.989e-04  -3.471e+00 
-3.378e4-00  -7.946e-04  -3.480ef00 
-3.341e+00  -1.297e-03  -3.443e+00 
-3.356e+00  1.524e-03  -3.457e+00 
-3.372e+00  1.48le-03  -3.471e+00 
-3.388e+00  -6.407e-04  -3.483e+00 
-3.341e+00  -1.365e-03  -3.434e+00 
-3.378ef00  -1.204e-03  -3.466ef00 
-3.385e+00  1.285e-03  -3.468e4-00 
-3.356e+00  1.706e-03  -3.439e+00 


TABLE  I.  Program  EAGLE  -  Flow  Solv 
NACA  0012  Airfoil  CXatput  (oDntinued/ 


52 

0.6698 

53 

0.6483 

54 

0.6266 

55 

0.6048 

56 

0.5828 

57 

0.5608 

58 

0.5387 

59 

0.5168 

60 

0.4949 

61 

0.4732 

62 

0.4517 

63 

0.4305 

64 

0.4097 

65 

0.3892 

66 

0.3692 

67 

0.3496 

68 

0.3305 

69 

0.3120 

70 

0.2940 

71 

0.2766 

72 

0.2598 

73 

0.2437 

74 

0.2281 

75 

0.2132 

76 

0.1989 

77 

0.1852 

78 

0.1722 

79 

0.1593 

80 

0.1480 

81 

0.1368 

82 

0.1262 

83 

0.1161 

84 

0.1066 

85 

0.0977 

8b 

0.0893 

87 

0.0813 

88 

0.0739 

89 

C.0669 

90 

0.0603 

91 

0.0542 

92 

0.0485 

93 

0.0432 

94 

0.0382 

95 

0.0336 

96 

0.0294 

97 

0.0254 

98 

0.0218 

99 

0.0185 

100 

0.0155 

101 

0.0127 

102 

0.0103 

103 

0.005 i 

104 

0.0062 

105 

0.0045 

106 

0.0031 

107 

0.0020 

108 

0.0012 

-0.1300 
-0.1486 
-0.1676 
-0.1871 
-0.2068 
-0.2264 
-0.2470 
-0.2683 
-0.2901 
-0.3121 
-0.3350 
-0.3589 
-0.3766 
-0.4030 
-0.4143 
-0.3353 
-0.5631 
-0.6568 
-0.6574 
-0.6523 
-0.6453 
-0.6357 
-0.6237 
-0.6094 
-0.5934 
-0.5758 
-0.5571 
-0.5372 
-0.5156 
-0.4920 
-0.4665 
-0.4425 
-0.4151 
-0.3877 
-0.3572 
-0.3289 
-0.2970 
-0.2617 
-0.2285 
-0.1907 
-0.1490 
-0.1112 
-0.0614 
-0.0159 
0.0383 
0.0942 
0.1567 
0.2226 
0.2964 
0.3741 
0.4632 
0.5573 
0.6579 
0.7591 
0.8617 
0.9511 
1 .0374 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (continued) 


109 

0.0006 

1 .0908 

no 

0.0002 

1.1088 

111 

0.0000 

1.1022 

112 

0.0000 

1 .0670 

113 

0.0002 

0.9998 

114 

0.000 

0.9000 

115 

0.0012 

0.7841 

116 

0.0020 

0.6487 

117 

0.0031 

0.5123 

1 18 

0.0045 

0.3758 

119 

0.0062 

0.2541 

120 

0.0081 

0.1394 

121 

0.0103 

0.0357 

122 

0.0127 

-0.0576 

123 

0.0155 

-0.1392 

124 

0.0185 

-0.2104 

125 

0.0218 

-0.2747 

126 

0.0254 

-0.3412 

127 

0.0294 

-0.3833 

128 

0.0336 

-0.4348 

129 

0.0382 

-0.4884 

130 

0.0432 

-0.5153 

131 

0.0485 

-0.5520 

132 

0.0542 

-0.5955 

133 

0.0603 

-0.6318 

134 

0.0669 

-0.6596 

135 

0.0739 

-0.6842 

136 

0.0813 

-0.7104 

137 

0.0893 

-0.7398 

138 

0.0977 

-0.7675 

139 

0.1066 

-0.7908 

140 

0.1161 

-0.8108 

141 

0.1262 

-0.8309 

142 

0.1368 

-0.8540 

143 

0.1480 

-0.8772 

144 

0.1598 

-0.8987 

145 

0.1722 

-0.9174 

146 

0.1852 

-0.9335 

147 

0.1989 

-0.9489 

148 

0.2132 

-0.9658 

149 

0.2281 

-0.9834 

150 

0.2437 

-1.0006 

151 

0.2598 

-1.0163 

152 

0.2766 

-1.0304 

153 

0.2940 

-1.0429 

154 

0.3120 

-1.0543 

155 

0.3305 

-1.0650 

156 

0.3496 

-1.0759 

157 

0.3692 

-1.0870 

158 

0.3892 

-1 .0979 

159 

0.4097 

-1.1079 

160 

0.4305 

-1.1169 

161 

0.4517 

-1.1246 

162 

0.4732 

-1.1310 

163 

0.4949 

-1.1364 

164 

0.5168 

-1.1410 

165 

0.5387 

-1.1452 
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TABLE  I.  Program  EAGLE  -  Flow  Solver 
NACA  0012  Airfoil  Output  (concluded) 


166 

0.5608 

-1.1486 

167 

0.5828 

-1.1557 

168 

0.6048 

-1.1557 

169 

0.6266 

-0.3054 

170 

0.6483 

0.0293 

171 

0.6698 

0.0177 

172 

0.6909 

0.0179 

173 

0.7117 

0.0202 

174 

0.7321 

0.0244 

175 

0.7521 

0.0307 

176 

0.7717 

0.0389 

177 

0.7907 

0.0486 

178 

0.8092 

0.0598 

179 

0.8271 

0.0720 

180 

0.8445 

0.0855 

181 

0.8612 

0.1003 

182 

0.8773 

0.1161 

183 

0.8929 

0.1334 

184 

0.9078 

0-1519 

185 

0.9220 

0.1731 

186 

0.9357 

0.1955 

187 

0.9487 

0.2215 

188 

0.9611 

0.2524 

189 

0.9729 

0.2924 

190 

0.9842 

0.3428 

191 

0.9948 

0.3753 

MACH  NO. 

PHIX 

PHIY  PHIZ 

CFL  LIFTAX 

0.8000 

0.00 

0.00  -1.25 

15.0 

2 

REFEI^ICE 

SURFACE  NAME 

AREA 

X 

Y 

Z 

1  NACA  0012 

Test 

1.00000 

0.25000 

0.00000  0.00000 

SURFACE  NAME 

FDRCIS 

MOMENTS 

1  NACA  0012 

Test 

FX  = 

0.009289 

MX  = 

-0.167376 

FY  = 

0.334751 

MY  = 

0.004645 

FZ  = 

0.000000 

FE  = 

0.031203 

COEFFICIENTS 

SURFACE  NAME 

FORCES 

MOMENTS 

1  NACA  0012  Test 

CX  = 

0.009289 

CMX  = 

-0.167376 

a  = 

0.334751 

CMY  = 

0.004645 

CZ  = 

0.000000 

CMZ  = 

0.031203 

SURFACE  NAME 

FORCES 

OOEFFTCIENTS 

o 

o 

to 

Test 

FL  = 

0.334469 

CL  = 

0.334469 

FT)  = 

0.016589 

CD  = 

0.016589 

FS  = 

0.000000 

CS  = 

0.000000 

TOIAL  CP  TIME  WAS  241.545929 
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2.  GENERIC  ETNNED  OCNFIGURATICN 

Ihis  configuration  is  a  generic  weapon  shape  with  fins 
consisting  of  a  tangent  ogive  fordbody  and  afterbody  with  a 
cylindrical  center  section  and  a  set  of  four  NACA  0008  airfoil 
fins.  This  airframe  has  been  specifically  designed  to  be  used 
for  CFD  code  validation.  Ihe  Air  Force  Armament  Laboratory  has 
expended  a  substantial  amount  of  resources  to  ccitpile  an 
extensive  surface  pressure,  force  and  mcment  data  base  for  this 
configuration  (Figure  13),  as  well  cis  the  same  configuration  in 
two  and  three  store  mutual  interference  combinations,  for  the 
validation  of  computational  aerodynamic  techniques  (Section  6.4). 

a.  Ebur  Block  Grid 

Ihe  computational  region  for  this  body  is  divided  into  four 
circumferential  blocks  (Figure  14a)  from  -45  to  45  degrees  (Block 
I),  45  to  135  degrees  (Block  II),  135  to  225  degrees  (Block  III), 
and  from  225  to  315  or  -45  degrees  (Block  IV).  The  entire  grid 
is  made  up  of  these  four  140  X  24  X  10  blocks  divided  equally 
between  the  four  fins  so  as  to  allow  the  solid  body  (impermeable 
wall)  boundary  conditions  to  be  implemented  on  a  block  boundary. 
As  before,  the  grid  lines  are  concentrated  down  along  the  solid 
body  (the  J=1  line)  and  along  the  leading  edge  of  the  fins  (the 
1=68  lines ) .  In  both  of  these  areas,  sharp  gradients  are  expected 
to  occur  in  the  flowfield  (Figure  14b). 

b.  Flow  Solver  Input 

The  input  for  the  generic  finned  configuration  will  specify 
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Figure  13.  Generic  Finned  Configiuration 
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Generic  Finned  Configuration  A  Block  Grid 
(Front  View) 


Figure  14b.  Generic  Finned  Configuration  4  Block  Grid 

(Side  View) 


91 


the  first  supersonic  test  case  (Mach  =  1.05)  at  a  fifteen  degree 
incidence  angle.  Here  we  have  iiipleinented  a  full  36CHdegree  grid 
so  as  to  allow  future  investigations  into  yaw  and  roll  angle 
effects. 

E$  FINPOT  CFL=5.0,  F91ACH=1.05,  JFRBQ=1,  SPLIT=2,  LIMIT=1, 
NSURF=4  ,  PC«DE3l=2,  NPR=1 ,  NIT=1000,  PRIWr=123, 
NGRAD=10.  hEPG=50  $ 

E$  POTATE  'IRANS=3,  PHIZ=  -15.0$ 

Body  w/sting 

E$  SURFACE  SREF=1,  XREF=  6.4  $ 

Fin  Set 

E$  SURFACE  SREF=2,  XREF=  0.25  $ 

Body  Alone 

E$  SURFACE  SREF=3,  XREF=  6.4  $ 

Generic  Store 

E$  SURFACE  SREF=4,  XREF=  6.4  $ 

As  before,  the  CEL  condition  is  set  first  (CEL=5.0)  along  with 
the  Mach  nunter  (FSMACE=1 .05)  and  the  nunfcer  of  updates  to  the 
flux  Jacobian  matrices  (JFE^=1 ).  For  this  specific  case,  the 
Jacobians  will  be  calculated  every  iteration.  A  form  of  Roe 
averaging  will  be  enployed  (SPLIT=2)  using  the  second-order 
(R0RDEE=2)  MINMOD  limiter  (LIMIT=1).  A  total  of  1000  iterations 
will  be  run  (NPR*NIT)  with  print  outs  occuring  once  every  1000 
iteraticxis  (NPR=1 ).  Ihe  print  format  will  be  x,  y,  and  z 
coordinates  versus  Cp  (pressure  coefficients),  as  well  as  all 
forces  and  moments  for  the  specified  4  surfaces  (NSURF=4).  Zero 
pressure  gradient  boundary  conditions  will  be  applied  for  the 
first  10  iterations  (NC3RAD=10)  with  a  total  of  50  zero  pressure 
gradient  boundary  conditions  (NZPG=50)  for  the  entire  run.  The 
rotations  for  this  case  will  be  fifteen  degrees,  and  therefore 
will  be  input  so  as  to  rotate  the  body  about  the  z  axis.  The 
surface  inputs  specify  four  different  sections  of  the  body  to 
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calculate  forces  and  moments  on  or  about.  As  is  always  the  case, 
the  first  portion  of  the  body  is  the  entire  configuration 
(SREF=1 )  entitled  "Body  w/sting".  The  second  section 
enooitpasses  the  entire  fin  set  (SREF=2),  entitled  "Fin  Set"  and 
is  specified  in  the  BCIN  input  cards  with  SURF=2.  The  third 
portion  is  ccmprised  of  the  ogive-cylinder-ogive  body  (SREF=3), 
called  "Body  Alone",  and  is  specified  in  the  BCIN  input  cards 
with  SURF=3.  And  finally,  the  fourth  specified  portion  is  the 
"Generic  Store"  which  is  ccirprised  of  just  the  body  and  fins  (not 
including  the  sting).  The  boundary  condition  inputs  are  as 
follows  for  this  four  circumferential  block  case: 

E$  BCIN  BCTYPE=1,  SURF=1 ,3,4,BLKA=1 ,  STARTA=  1,  1,  1, 

ENDA=  116,  1,  10  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,3,4,BLKA=2,  STAFTA=  1,  1,  1, 

ENnA=  116,  1,  10  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,3,4,BLKA=3,  STARTA=  1,  1,  1, 

ENDA=  116,  1,  10  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,3,4,BLKA=4,  STARTA=  1,  1,  1, 

ENDA=  116,  1,  10  $ 

BCIN  BCTYPE=1,  SURF=1,  BLKA=1,  STARTA=1 17, 1 ,  1, 

ENnA=  140,1  ,  10$ 

E$  BCIN  BCTYPE:=1,  SURF=1,  BLKA=2,  STAPTA=117,1,  1, 

ENnA=  140,1  ,  10$ 

E$  BCIN  BCrYPE=1,  SURF=1,  BLKA=3,  STARTA=1 17, 1 ,  1, 

ENnA=  140,1  ,  10$ 

E$  BCIN  BCTYPE=1,  SURF=1,  BLKA=4,  STAKI!A=1 17, 1 ,  1, 

ENnA=  140,1  ,  10$ 

E$  BCIN  BCTYPE=1,  SURF=1 ,2,4,BLKA=1 ,  STARTA=68,  1,1, 

ENDA=  113,16,  1  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,2,4,BLKA=1 ,  STARTA=68,  1,10, 

ENDA=  113,16,  10  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,2,4,BLKA=2,  STAREA=68,  1,1, 

ENnA=  113,16,  1  $ 

E$  BCIN  BCTTOE=1,  SURF=1 ,2,4,BLKA=2,  START!A=68,  1,10, 

ENDA=  113,16,  10  $ 

E$  BCIN  BCTyPE=1,  SURF=1 ,2,4,BLKA=3,  STAKIA=68,  1,1, 

ENDA=  113,16,  1  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,2,4,BLKA=3,  STAKEA=68,  1,10, 

ENDA=  113,16,  10  $ 

E$  BCIN  BCTn>E=1,  SURF=1 ,2,4,BLKA=4,  STARTA=68,  1,1, 

ENDA=  113,16,  1  $ 

E$  BCIN  BCTYPE=1,  SURF=1 ,2,4,BLKA=4,  STARTA=68,  1,10, 

ENDA=  113,16,  10  $ 

E$  BCIN  BCTYPE=2,  BLKA=1,  STARTA=  1,24,  1, 

ENDA=  140,24,  10  $ 

E$  BCIN  BCTYPE=2,  BLKA=1,  STARTA= 1 40 , 1 ,  1, 

ENDA=  140,24,  10  $ 
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E$  BC3N  BCTYPE=2, 

BLKA=2,  STARTA=  1,24,  1, 
ENDA=  140,24,  10  $ 

E$  BCm  BCTYPE=2, 

BLKA=2,  STARTA=140, 1,  1, 
ENDA=  140,24,  10  $ 

E$  BCIN  BCTYPE=2, 

BLKA=3,  STARIA=  1,24,  1, 
ENDA=  140,24,  10  $ 

E$  BCIN  BCTYPE=2, 

BLKA=3,  STARI!A=140,1 ,  1, 
ENnA=  140,24,  10  $ 

E$  BCIN  BCT7PE=2, 

BLKA=4,  STARrA=  1,24,  1, 
ENnA=  140,24,  10  $ 

E$  BCIN  BCIYPE=2, 

BLKA=4,  STARTA=140,1,  1, 
ENDA=  140,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=1,  STARTA=  1,  1,  1, 
ENDA=  68,24,  1  , 

BLKB=4,  STARrB=  1,  1,10, 
ENrB=  68,24,  10  $ 

E$  BCIN  BCTYPE:=4, 

BLKA=2,  STAm=  1,  1,  1, 
ENDA=  68,24,  1  , 

BLKB=1,  STARTB=  1,  1,10, 
ENDB=  68,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=3,  STARTA=  1,  1,  1, 
ENDA=  68,24,  1  , 

BLKB=2,  STAKIB=  1,  1,10, 
ENDB=  68,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=4,  STARIA=  1,  1,  1, 
ENnA=  68,24,  1  , 

BLKB=3,  STARrB=  1,  1,10, 
ENDB=  68,24,  10  $ 

E$  BCIN  BdYPE:=4, 

BLKA=1,  STARrA=68,16,  1, 
ENnA=  113,24,  1  , 

BLK&=4,  STARIB=68,16,10, 
ENDB=  113,24,  10  $ 

E$  BCIN  BdYPE=4, 

BLKA=2,  STARIA=68,16,  1, 
ENDA=  113,24,  1  , 

BLKB=1,  STARTB=68,16,10, 
ENDB=  113,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=3,  STARIA=68,16,  1, 
ENDA=  113,24,  1  , 

BLKB=2,  STARTB=68,16,10, 
ENDB=  113,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=4,  STARrA=68,16,  1, 
ENDA=  113,24,  1  , 

BLKB=3,  STARrB=68,16,10, 
ENDB=  113,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=1,  STARTA=113,1,  1, 
ENDA=  140,24,  1  , 

BLKB=4,  STARTB=113,1,10, 
ENDB=  140,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=2,  STARTA=113,1,  1, 
ENDA=  140,24,  1  , 

BLKB=1,  STARTB=1 13,1,10, 
ENDB=  140,24,  10  $ 

E$  BCIN  BCTyPE=4, 

BLKA=3,  STARTA=1 13,1 ,  1, 
ENDA=  140,24,  1  , 

BLKB=2,  STARTB=1 13,1,10, 
ENDB=  140,24,  10  $ 

E$  BCIN  BCTYPE=4, 

BLKA=4,  STARTA=I 13,1,  1, 

a'JDA=  140,24,  1  , 

BLKB=3,  STARTB=113,1,10, 
ENDB=  140,24,  10  $ 


E$  BCIN  BCTYPE=5, 


E$  BCIN  BCTYPE=5, 


E$  BCIN  BCIYPE=0  $ 


BLKA=1,  STARIA=  1,  1,  1, 
ENDf.=  1,24,  10  , 
BLKB=3,  STARTB=  1,  1,  1, 
ENDB=  1,24,  10  $ 
BLKA=2,  STARIA=  1,  1,  1, 
ENDA=  1,24,  10  , 
BLKB=4,  STARIB=  1,  1,  1, 
ENDB=  1,24,  10  $ 


■nie  boundary  condition  input  for  this  four  block  case  is  set  up 
for  ease  of  understanding  since  the  bc's  are  essentially 
syrtinetric.  The  first  four  inputs  are  for  the  solid  body 
(BCTYPE=1  )  ccxidition  along  the  J=1  surface  in  blocks  1-4  for  the 
body  alone.  These  inputs  correspond  to  the  SREF=1,  3  and  4  cards 
in  the  SURFACE  input.  Ihe  second  four  cards  specify  the  sting 
portion  of  the  configuration  and  correspond  to  the  SREF=1  ir^xit 
in  the  SURFTiCE  section.  The  next  eight  input  cards  are  for  the 
upper  and  lower  surfaces  of  the  fins  in  all  four  blocks  (K=1  and 
10).  These  are  also  inpermeable  wall  conditions  (BCrYPE=1 )  and 
correspcnd  to  the  SURFACE  input  SREF=1,  2  and  4.  Far field 
conditions  (BCTYPE=2)  are  then  inplemented  for  all  four  blocks  at 
the  outer  (J=24)  and  back  (1=140)  boundaries.  Block-to-block 
conditions  are  also  inplemented  (BCIYPE=4),  at  the  interface 
between  the  blocks  at  the  -45,  45,  135,  and  225  degree  positions 
(K=1  and  10).  The  final  two  input  cards  involve  the  stagnation 
line  out  the  front  of  the  configuration,  along  which  a  synmetry 
conditic*!  must  be  inplemented  (BCTYPE=5).  This  aids  in  resolving 
the  flux  through  what  is  effectively  the  zero  area  cells  on  the 
1=1  surface.  And  finally,  the  input  stream  is  concluded  with  a 
BCIYPE=0  card  terminating  all  input  for  this  run.  The 
conputational  grid  is  shown  in  Figure  15  to  give  the  proper 
ijTplementaticxi  of  the  boundary  condition  types  (BCTYPEs). 
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c.  Flow  Solver  Output 

TTie  output  for  the  Generic  Store  configuration  shows  the 
entire  run  took  10077.12623  CP  seconds  for  1000  iterations.  Ihe 
density  residuals  converged  almost  three  orders  of  magnitudes  as 
did  the  L2  Norm  depending  on  block  number.  The  nurrber  of 
supersonic  points  appears  to  be  "frozen"  for  the  last  20 
iteraticais.  This  seems  to  be  a  fairly  well  converged  solution 
and  should  yield  excellent  engineering  approximations  to  the  flow 
at  these  conditions. 
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SECTION  VI 


RESULTS 

This  sec±ion  surmvarizes  the  results  obtained  with  the  EAGLE- 
Flow  Solver  over  the  past  few  nonths  of  development  and 
validation.  The  code  has  been  exercised  thoroughly  on 
configurations  ranging  from  the  NACA  0012  airfoil  (using  grids 
created  by  the  EAGLE  -  Numerical  Grid  Generation  System  as  well 

as  the  GRAPE  code  developed  by  Sorenson  to  the  three  store 
mutual  interference  configuration  (Section  6.4).  An  entire  range 
of  Mach  numbers  (from  0.40  to  6.81 )  and  angles  of  attack  (from  0 
to  20  degrees)  have  been  run  to  test  the  code's  ability  to  solve 
for  difficult  flow  conditions.  This  section  provides  just  a 
sanpling  of  these  results. 

1.  NACA  0012  AIRPDIL 

The  NACA  0012  airfoil  has  been  used  as  a  test  configuration 
due  to  its  vast  amount  of  documentation  by  other  code  development 
efforts  and  because  it  is  so  well  understood  from  a  theoretical 
stanc^ioint.  The  test  case  run  in  Section  5.1  yields  excellent 
results  for  the  flow  conditions  spiecified  (Mach=  0.80,  A0A=  1.25 
degrees).  The  code  captures  both  the  upper  surface  shock  at  the 
55%  chord  location,  as  well  as  the  weak  lower  surface  shock 
present  at  the  35%  chord  location  (Figure  16).  No  experimental 
data  is  shown;  however,  the  solution  conpares  favorably  to  other 
codes results. 
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Figure  16.  NACA  0012  Airfoil  -  Flow  Solver  Output 
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2.  GENERIC  FINNED  OONFIGURATI»I 

The  Generic  Finned  Configuration  (Ogive-Cylinder-Ogive  with 
fins)  has  been  exercised  for  several  transonic  Mach  numbers  (0.80 
through  1.20)  and  angles  of  attadc.  This  configuration  has  been 
the  "bread  and  butter"  test  case  for  the  Armament  Laboratory  for 
freestream,  and  interference  flow,  analysis  in  the  transonic  Mach 
1 9  7.0 

range. The  following  results  are  provided  at 

Mach  1.05  at  two  angles  of  attack  to  show  the  codes  strengths  and 

wealcnesses  cis  flow  conditions  become  increasingly  viscous  in 

nature. 

Figure  17  shows  the  inviscid  solutions  obtained  for  the 
sanple  case  presented  in  Section  5.2  (Mach=  1.05  and  A0A=  15 

degrees)  versus  experimental  data^^.  In  Figure  17a,  the  body 
pressure  distribution  is  shown  at  an  orientation  angle  of  0 
degrees  (or  twelve  o'clock-  the  leeward  side).  Hie  solver  does 
an  excellent  job  of  following  the  flow  as  it  accelerates  over  the 
ogive  nose  and  e>q)ands  near  the  shoulder  region  of  the  body.  The 
inviscid  approximation  does  begin  to  fail  at  about  80%  as  the 
flow  separates  over  the  ogive  boattail  section.  Figure  17b 
shows  the  body  pressures  for  the  windwcud  side  (orientation  of 
180  degrees  or  six  o'clock  position).  The  flow  solver  again 
captures  the  flow  as  it  expands  over  the  nose  and  does  a  good  job 
predicting  the  flow  in  the  tail  section  (as  the  flow  stays 
attached  on  the  windward  side).  In  Figure  17c  fin  pressures  are 
shown  versus  experimental  data  at  a  fin  orientation  of  45  degrees 
and  a  span  location  of  20%  on  the  leeward  side  of  the  fin.  As 
expected,  the  inviscid  solver  doesn't  capture  the  physics  of  the 
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Figure  17a.  Generic  Finned  Configuration  -  Flow  Solver  Output 
Body  Pressure  Distribution  (  oi  =  15  Degrees) 
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Figure  17c.  Generic  Finned  Configuration  -  Flow  Solver  Output 
Fin  Pressure  Distribution  (  “  =  15  Degrees) 
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Figure  17d.  Generic  Finned  Configuration  -  Flow  Solver  Output 
Fin  Pressure  Distribution  (  “  =  15  Degrees) 


flcwfieM  in  an  area  where  sepenrated  flow  predGminates  the 
flowfieM  (leeward  side  of  the  body  on  the  leeward  side  of  the 
fin).  Figure  17d  shows  the  windward  side  of  the  same  fin  (fin 
orientation  of  45  degrees)  at  a  span  location  near  80%.  The 
solver  does  a  better  jc±)  of  capturing  the  expansion  over  the 
leading  edge  amd,  because  of  a  lack  of  data  points,  there  is  no 
judging  the  trailing  edge. 

Figure  18  gives  results  detained  by  the  flow  solver  at  Mach 
=  1.05  and  an  angle  of  attack  of  6  degrees.  The  code  does  an 
excellent  jcto  of  capturing  the  expansion  over  the  nose  on  to  the 
shoulder  region  of  the  body  at  the  leeward  (0  deg. ),  windward 
(180  deg.)  and  the  side  meridian  (90  deg.)  positions  (Figure  18a). 
Figure  18b,c  shows  results  dotained  for  the  windward  fins  (at  135 
deg. )  on  the  leeward  side  at  four  fin  span  locations  and  for  the 
same  fin  ( 1 35  deg. )  on  the  windward  side  at  four  fin  span 
locations.  In  both  cases  the  code  does  a  good  job  of  predicting 
the  surface  pressures  along  the  chord  length  of  the  fin. 

3.  BLONT,  FINNED  BODY  OF  REVOLOTICN  WITH  (»NARDS 

This  configuration  has  been  used  as  a  test  case  for  the 
EAGLE  -  Numerical  Grid  Generation  Systan,  as  well  as  for  the 
EAGLE  -  Flow  Solver,  due  to  its  difficult  canard  set  and  multi¬ 
block  layout  shown  in  Reference  21 .  Several  studies  have  been 
acoonplished  on  this  body  due  to  its  ability  to  create  a  multiple 

shock  structure.  '  Tne  results  shown  here  in 

Figure  19  are  for  a  Mach  number  of  0.95  and  angle  of  attack  of 

zero  degrees  for  both  the  EAGLE  -  Flow  Solver  and  a  previous 
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Figure  18a.  Generic  Finned  Configuration  -  Flow  Solver  Output 
Body  Pressure  Distribution  ( a  =  6  Degrees) 
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Figure  18b.  Generic  Finned  Configuration  -  Flow  Solver  Output 
Fin  Pressure  Distribution  ( “  =6  Degrees) 


Figure  19b.  Blunt,  Finned  Body  of  Revolution  with  High  Taper 
Ratio  Canards 


euqjlicit  algorithm  that  was  a  candidate  for  inclusion  in  the 

code^*^.  Both  codes  do  an  excellent  job  of  following 
the  e)q>ansion  over  the  blunt  nose  and  past  the  canards  at  the  15% 
-  20%  locations.  The  solver  predicts  the  ccnplex  flow  pattern 
throughout  the  rest  of  the  field,  but  fails  at  the  boattail 
section  as  the  flow  separates  off  the  back  end  (Figure  19a).  In 
Figure  19b,c  canard  pressures  are  shown  for  two  span  locations 
(30%  and  60%).  In  both  cases  the  codes  do  well  in  following  the 
pressure  distribution  over  the  canard  shape. 

4.  MLTIUAL  INTKKFEKENCE  COIFIGURATTCNS 

Tb  illustrate  the  capability  for  predicting  interference 
flew  conditicais,  the  flow  solver  was  run  for  a  set  of 
configurations  shown  in  Figure  20.  Ihe  geometry  is  that  shown  in 
Figure  14  and  a  discussion  of  the  single  and  triple  store  grids 
can  be  found  in  Volume  I. 

Solutions  were  obtained  for  the  three  store  configurations  at 
two  Mach  nuirbers  (0.80  and  0.95)  at  zero  and  ten  degrees  angle  of 
attack.  For  the  one,  two  and  three  store  cases.  Figure  21  plots 
the  surface  pressure  distribution  on  the  inboard  side  of  the  body 
(265  deg.  location).  As  expected,  the  interference  increases 
as  the  second  and  third  bodies  are  included.  The  aerodynamic 
interference  manifests  itself  in  stronger  expansions  just  after 
the  shoulder  region  that  increases  in  size  downstream  as  the 
other  bodies  are  added.  The  flow  solutions  agree  quite  well  with 
the  ejqjerlmental  data;  Ixjwever,  the  code  over  predicts  the 
ejgansions  at  the  shocks  as  well  as  the  shock  locations. 
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Figure  20.  fftiltiple  Body  Configurations 
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To  investigate  their  effect,  calculations  were  c±itained  for 
both  the  finned  and  unfinned  configurations.  Figure  22  shows 
the  effects  of  the  fins  on  the  double  (or  two  store) 
configuration  at  Mach=  0.95  and  zero  degrees  angle  of  attack. 

The  effect  seen  on  the  body  in  the  interference  flow  region  is  a 
reductican  of  the  expansion  just  aft  of  the  nose  for  the  finned 
body,  but  an  increased  expansion  in  the  fin  area.  The  fins  also 
increase  blockage  in  the  flow,  moving  the  fordxdy  shock  forward 
3%. 

To  examine  the  interference  effects  at  angle  of  attack,  the 
triple  store  configuration  with  fins  was  run  for  both  zero  and 
ten  degrees  angle  of  attack  at  Mach=  0.80.  The  surface  pressure 
distribution  for  these  cases  is  shown  in  Figure  23.  As  the 
configuration  increases  in  incidence  angle,  the  fordbody 
expansion  decreases  on  the  leeward  side  of  the  lower  body.  The 
upper  two  stores  in  the  three  store  configuration  act  to  redirect 
the  flow  over  the  lower  body  and  effectively  reduce  the  angle  of 
attack  seen  by  the  lower  body.  This  accounts  for  the  weaker 
forebody  expansion  seen  at  the  ten  degree  incidence  angle. 

The  interference  effect  of  the  fins  was  investigated  at  Mach= 
0.80  and  zero  degree  incidence.  Figure  24  illustrates  the 
pressure  distributions  on  the  fins  for  the  three  store 
configuration  at  the  positions  of  minimum  and  maximum 
inter  femce.  The  pressure  distribution  on  the  minimum 
interference  fin  is  nearly  flat  while  a  significant  expansion  and 
shock  occurs  on  the  maximum  interference  fin.  This  is  expected 
since  the  flow  accelerates  through  the  region  between  the  three 
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bodies  lowering  the  pressure.  Ihe  flow  solver  does  a  fair  to 
good  job  on  the  first  50%  of  the  chord  for  both  fins.  Ihe  shock 
on  the  inboard  fin  is  missed  substantially  both  in  strength  and 
location. 
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SECTIOK  VII 


SUMMARY 

As  presented  here,  the  Program  EAGLE  -  Flow  Solver  can  be 
enployed  cis  an  extremely  effective  aid  in  the  analysis  of 
freestream  and  interference  flow  aerodynamic  characteristics  of 
advanced,  arbitrarily  shaped  airframe  configurations. 

The  code  appears  to  be  most  effective  in  the  high  subsonic 
to  low  supersonic  Mach  range  (0.75  to  2.0),  but  has  shown 
excellent  results  outside  of  this  range.  As  is  evident  from  the 
inviscid  nature  of  the  code,  the  solver  has  difficulty  resolving 
viscous  dominated  regions  (separated  flow)  and;  therefore  is 
limited  to  relatively  small  incidence  angles.  However,  the  flow 
solver  has  shown  excellent  results  on  body  pressure  distributions 
at  incidence  angles  as  high  as  15  degrees. 

Good  engineering  solutions  can  also  be  ctotained  for  cottplex 
geometries  in  interference  flowfields  as  is  evident  from  Section 
6.4.  The  solutions  cttained  show  that  an  inviscid  approximation 
is  quite  adequate  for  most  aspects  of  the  fluid  dynamic 
phenomena;  however,  the  solution  does  break  down  in  highly 
viscous  regions  such  as  the  p»sitions  of  maximum  interference 
where  the  boundary  layer  is  "sucked"  back  into  the  channel 

between  the  stores. 

In  general,  the  EAGLE  -  Flow  Solver  is  an  excellent 
analytical  tool  for  the  practicing  engineer  in  r,Ged  of  a 
cotpliment  to  eaq^erimental  analysis.  The  EACT.E  code  should  yield 
extremely  accurate  predictions  of  the  aerodynamic  characteristics 


120 


of  any  arbitrarily  shaped  advanced  airframe  configuration  in 
freestream  or  interference  flowfields.  Tbe  recoirmended  range  of 
inplementation  is  for  subsonic  to  low  hypersonic  Mach  numbers 
(where  there  are  no  real  gas  effects  present)  [  0.40  <  Mach  <  8.0 
at  relatively  low  angles  of  incictence  (where  the  s^>arated 
region  is  relatively  small)  [  0  to  20  degrees  ]. 
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